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Results of Current Research 

1) The verification of LAWS beam alignment in space can be achieved by a measurement of heterodyne 
efficiency using the surface return. Due to the statistical nature of the return signal, about 40 shots will be 
required for a good test This verification method does not require knowledge of the backscatter 
coefficient of the surface, the power transmitted, the atmospheric attenuation, or the detector gain. The 
crucial element is a direct detection signal that can be identified for each surface return. This should be 
satisfied for LAWS but will not be satisfied for descoped LAWS. These results were presented at the 
Coherent Laser Radar Meeting in Paris. The behavior of estimates for heterodyne efficiency in the near- 
field (ground based testing of LAWS) were also determined. For the actual LAWS optical geometry, the 
required sensitivity is not available. A focused geometry is required to test the space based system. 

2) The performance of algorithms for velocity estimation can be described with two basic parameters: the 
number of coherently detected photo-electrons per estimate and the number of independent signal sam- 
ples per estimate. For low signal levels, the fraction of bad estimates is required to describe performance. 

The fraction of bad estimates is mainly a function of the number of effective photo-electrons per estimate 
and has little dependence on the number of independent samples per estimate; is independent of the esti- 
mator used; and has a slight dependence on the velocity space searched. 

3) The average error of spectral domain velocity estimation algorithms are bounded by a new periodo- 
gram Cramer-Rao Bound. Comparison of the periodogram CRB with the exact CRB indicates a factor of 
two improvement in velocity accuracy is possible using non-spectral domain estimators. This improve- 
ment has been demonstrated with a maximum-likelihood estimator. 

4) The comparison of velocity estimation algorithms for 2 and 10 micron coherent lidar was performed by 
assuming all the system design parameters are fixed (range to target, telescope aperture, detector quantum 
efficiency, heterodyne efficiency, perfect beam alignment, optical element efficiency, laser pulse energy, 
atmospheric attenuation, range resolution, observation time per estimate, velocity search space) and the 
signal s tatis tics are dominated by a 1 m/s rms wind fluctuation over the range gate. Then the quality of 
the data product is only a function of the wavelength dependence of the backscatter coefficient (3. For 
pocjf 1 the number of photons collected per estimate is the same and a 10 micron system has better data 
quality because the fraction of bad estimates at low signal levels is less than for the 2 micron lidar. For 
(JocXr 2 the signal to noise ratio is the same and the 2 micron lidar has better data quality. The two systems 
have similar data quality when |3«X -1 ' 3 . These results agree with the GE and TRW reports. 

5) The beam alignment requirements for 2 micron are much more severe than for a 10 micron lidar. The 
effects of the random backscattered field on estimating the alignment error is a major problem for space 
based lidar operation, especially if the heterodyne efficiency cannot be estimated. This is a sever handi- 
cap for 2 micron coherent lidars because the beam alignment must be changed by 5-6 microradians per 
minute due to variations in orbit height and 6 microradians of misalignment corresponds to more than 20 
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6) For LAWS, the biggest science payoff would result from a short transmitted pulse, on the order of 0.5 
microseconds instead of 3 microseconds. The advantages are 

a) 50-100 meter range resolution in the boundary layer, thus providing data of great interest to atmos- 
pheric scientists. 

b) More efficient laser design for both C02 and solid state. 

c) QUICK LAWS will only measure the boundary layer. This design is ideal for QUICK LAWS and can 
be scaled to full LAWS with minimal changes. 

d) For conditions of higher backscatter, the 50 m range resolution will permit useful measurements of 
wind statistics that are essential for optimal design of velocity estimators. 

e) A short pulse will permit better velocity estimation algorithms because the statistics of the signal will 
be determined by the transmitted pulse (which is known) instead of the velocity fluctuations over the 
sensing volume of the pulse (which are unknown and difficult to estimate). 

0 Better measurements near clouds and the surface will be possible. 

The disadvantage are: 

a) The velocity accuracy in the regions of low backscatter with 1 km height resolution will be about 1.5- 
2 m/s. 

b) A design of a high-energy short pulse C0 2 laser would be difficult to verify with the bread-board 
laser. 

7) The numerically errors for simulation of laser propagation in the atmosphere have been determined as 
a joint project with the University of California, San Diego. Useful scaling laws were obtained for Kol- 
mogorov atmospheric refractive turbulence and a atmospheric refractive turbulence characterized with an 
inner scale. This permits verification of the simulation procedure which is essential for the evaluation of 
the effects of refractive turbulence on coherent Doppler lidar systems. 

8) The analysis of 2 micron Doppler lidar data from Coherent Technologies, Inc. (CTI) has demonstrated 
many of the advantages of doppler lidar measurements of boundary layer winds. The estimation error for 
a single shot with 50 meter range resolution has been estimated at 0.65 m/s., which is within 40% of the 
ideal performance of the Cramer-Rao Bound. The effects of wind shear and wind turbulence over the 
pulse volume are probably the dominant source of the reduced performance. 

9) The effects of wind shear and wind turbulence on the statistical description of doppler lidar data has 
been derived and calculated. This is the first analysis for the single-shot statistics which are required for 
most Doppler lidar applications. The data is non-stationary and new algorithms are required to take full 
advantage of real Doppler lidar data. The effects of wind tuibulence over the range resolution cell are 
noticeable for typical LAWS parameters. The effects are more pronounced in the planetary boundary 
layer. The effects of wind shear and wind turbulence are more pronounced for 2 microns than 10 
microns. 
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Future Work 

1) Determine the ability of velocity estimation algorithm’s to identify which estimates are good when the 
signal is low. 

2) Include velocity fluctuations and wind shear in velocity estimators for LAWS. 

3) Include the estimation of signal parameters like signal power and spectral width in velocity estimators. 

4) Characterize the performance of velocity estimators for better OSSIE’s. 

5) Investigate beam alignment techniques for LAWS including the random nature of the received signals. 

6) Collect information about surface return statistics and their effects on beam alignment. 

7) Collect information about the wind statistics and their effects on velocity estimation. 
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Abstract- The performance of mean-frequency estimators for Doppler radar/lidar measurements of winds 
is presented in terms of two basic parameters: 0 the ratio of the signal energy per estimate to the spectral 
noise level and £2 which is proportional to the number of independent samples per estimate. For fixed <1> 
and £2, the the Cramer-Rao Bound (theoretical best performance) for unbiased estimators of mean- 
frequency (normalized by the spectral width of the signal), signal power, and spectral width, are essen- 
tially independent of the number of data samples M . For large 0 , the estimators of mean-frequency are 
unbiased and the performance is independent of M. The spectral domain estimators and covariance based 
estimators are bounded by the approximate periodogram Cramer-Rao Bound. The standard deviation of 
the Maximum-Likelihood estimator approaches the exact CRB, which can be more than a factor of two 
better than the performance of the spectral domain estimators or covariance based estimators for typical 
£2. For small 0 the estimators are biased due to the effects of the uncorrelated noise (white noise) which 
results in uniformly distributed "bad” estimates. The fraction of bad estimates is a function of 0 and M 
with weak dependence on the parameter £2. Simple empirical models describe the standard deviation of 
the good estimates and the fraction of bad estimates. For Doppler lidar and for large 0 , better perfor- 
mance is obtained by using many low energy pulses instead of one pulse with the same total energy. For 
small O, the converse is true. 


1. Introduction 

Measurement of wind fields using Doppler radar have been routinely performed for many years. 
The design and operation of Doppler radar is reviewed by Doviak and Zmic (1984, 1993). The optical 
counterpart to Doppler radar, coherent Doppler lidar, has become more important with the development 
of efficient C0 2 [Bilbro et al„ 1986; Petheram et al„ 1989; Post and Cupp, 1990; Pearson and Rye, 1992; 
Gal-chen et al„ 1992] and solid state lidars [Kavaya et al, 1992; Henderson et al., 1992, 1993]. The 
scatterers for Doppler radar are refractive index fluctuations, hydrometeors, and insects. The scatterers for 
Dopper lidar are atmospheric aerosol particles. Doppler radar data is generated using a sequence of many 
pulses, usually separated by 0.5-5 msec. Doppler lidar data is produced with a single pulse which permits 
many more estimates in a given time. The statistical description of the data is the same (ignoring ground 
clutter for radar) because both signals are produced by the super-position of many randomly phased scat- 
tered fields. Both of these instruments estimate the radial velocity of the scatterers from the Doppler fre- 
quency or mean-frequency of the signal by using various estimators. 

Most mean-frequency estimators are either spectral domain or covariance estimators, [Levin. 1965; 
Zmic, 1979; Mahapatra and Zmic, 1983; Miller and Rochwarger, 1972; Doviak and Zmic, 1984, 1993; 
May and Strauch, 1989; May et al., 1989; Rye and Hardesty, 1993a,b,c; Sirmans, and Bumgarner, 1975], 
i.e., they estimate the mean-frequency of the weather signal using estimates of the spectrum or estimates 
of the covariance [see also Kay and Marple, 1981; Marple, 1987]. 

The ideal performance of an unbiased estimator is given by the Cramer-Rao Bound (CRB) [Hel- 
strom, 1968; Van Trees, 1968; Scharf, 1991; Frehlich, 1993a]. If an unbiased estimator approaches the 
CRB, it is a Maximum-Likelihood (ML) estimator. Approximate calculations of the CRB have been dis- 
cussed by Whittle (1953), Levin (1965), Zmic (1979),Rye and Hardesty, (1993a,b). These approxima- 
tions are based on the spectrum of the signal, and therefore require a sufficiently long data sequence for 
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the spectrum to be well defined [see Marple, (1987) for a discussion of the effects of finite data and 
discrete sampling]. The calculation of the exact CRB has been discussed by Frehlich (1993a) who intro- 
duced an approximate CRB based on the average periodogram instead of the spectrum. We will compare 
the performance of various mean-frequency estimators to the exact CRB and the periodogram CRB. The 
CRB and the performance of mean-frequency estimators depend on the parameters of the data. The sys- 
tem parameters are chosen to emphasize the most important physical mechanisms of the problem. This 
reduces the parameter space to two basic variables: <t> the ratio of the signal energy per estimate to the 
spectral level of the noise and Q which is proportional to the number of independent samples per esti- 
mate. (For Doppler lidar <I> is also the number of effective photo-electrons coherently detected per esti- 
mate). For fixed O and fl, the number of data points M per estimate has little effect on performance, 
especially for large <D. This simplifies the design and analysis of Doppler radar/lidar systems. 

The performance of mean-frequency estimators has been described by it’s standard deviation [Zmic, 
1979; Sirmans, and Bumgarner, 1975; Mahapatra and Zmic, 1983; Doviak and Zmic, 1984, 1993; May 
and Strauch, 1989; May et al., 1989; Rye and Hardesty, 1993a,b] and the number of estimates that fall 
within a given region [Anderson, 1991; Rye and Hardesty, 1993a]. We describe the performance of 
mean-frequency estimators using an empirical model for the Probability Density Function (PDF) of the 
estimates. This will permit a meaningful comparison to the CRB when the estimates are biased due the 
effects of the uncorrelated noise that produces uniformly distributed "bad" estimates at low signal levels. 

For Doppler lidar. Rye and Hardesty, (1993a) considered the question: are many pulses of low 
energy better than one pulse of the same total energy? This question will be addressed for the unbiased 
and biased regimes. 

2. Basic System Parameters 

Doppler radars typically employ a complex receiver (mixers that convert the radio signals to com- 
plex data) to measure positive and negative velocities [Doviak and Zmic, 1984, 1993]. We consider com- 
plex data throughout this paper with specific references to the case of real data. Doppler radar or Doppler 
lidar data is well approximated as a zero-mean complex Gaussian random vector z with elements z k 
which satisfy <z k z, >=0, where <•> denotes ensemble average [Doviak and Zmic, 1984, 1993; Helstrom, 
1968]. The PDF of the data is Gauusian and defined by the covariance matrix R with elements 

=<z*z/* > ( 1 ) 

We assume the signal model 

z* = s k cxp(2KikfT s }¥n k (2) 

where f[Hz] is the mean frequency, T s [sec ] is the sampling interval, the random signal s k is indepen- 
dent of the uncorrelated noise n k where 

<n k nj>=0, <n k nJ>=N5 k .j , ( 3 ) 

N=<\n | 2 > is the average noise power, | n | denotes the modulus of the complex variable n , and 8 is the 
Kronecker delta symbol. 

For Doppler lidar, the mean frequency / is defined by the radial component of the velocity of the 
scatterers in the center of the sensing volume of a given range gate [Frehlich, 1993b]. The mean- 
frequency is well defined for a single shot for any random velocity field and random collection of scatter- 
ers. For Doppler radar, the mean-frequency is related to the spectrum of the wind field and sensing 
volume of the pulse [Doviak and Zmic, 1984, 1993]. For both Doppler radar and lidar, the signal model 
of Eq. (2) is valid. An estimate of the mean-frequency is produced with M data samples which defines 
the total observation time per estimate as MT S [sec ]. 

The average noise power is set to unity to simplify the results. The data is assumed stationary, i.e., 
Rki^k-i and the covariance reduces to 

R k = SNR pi exp[2 nifkT s ]+S* (4) 

where SNR =S IN is the Signal-to-Noise Ratio, S=< \ s | 2 > is the average signal power, and p* is the nor- 
malized covariance of the signal. The performance of mean frequency estimators will be presented for 
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the Gaussian covariance model of Zmic (1979) 

p* = sxvi-lT&iwkTs) 1 ] (5) 

where w [Hz ] is the spectral width of the signal. With this model, the data is fully characterized by the 
system parameters (f *S ,w) as well as the experimental parameters M and 7$ . For Doppler lidar, this sig- 
nal model is exact [Frehlich, 1993b] if there are no wind fluctuations over the sensing volume and a 
Gaussian laser pulse with temporal power profile Pi (f )[W ] given by 

F L (r)=exp(-r 2 /a 2 ) (6) 

is transmitted through a telescope aperture where a[sec ] is 

0=\/('$KW ) . ( 7 ) 

The Full- Width at Half-Maximum (FWHM) At [sec ] of the transmitter pulse is 

Ar=2Vln2o . (8) 

The 1/e Full Width of the signal covariance At [sec], which defines the signal correlation time (see Fig. 
1). is 

At=-^— = 4<t . (9) 

KW 

For a Doppler Lidar operating under ideal conditions (the detector noise is dominated by the local 
oscillator shot noise, the shot noise is determined by the Poisson statistics of the detected photons, and 
there is no atmospheric refractive turbulence, Frehlich and Kavaya, 1991) the SNR is given by 

SNR(t)=^~-r\ H (0 ( 10 > 

where P D (r) is the effective direct detection power measured by the detector, r| w (r) is the heterodyne 
efficiency, h =6.626 x 10 -3 Vs ] is Planck’s constant, v[//z] is the optical frequency, and B [Hz ] is the 
detector bandwidth. For ideal operation 

„ 00 

SNR ( t ) = t \\P L (f-2 R lc )[K (R )] 2 p(tf )C(R)dR (11) 

h\B J 0 

where r\ Q [electrons /photon ] is the detector quantum efficiency, c [m s _1 ] is the speed of light in a homo- 
geneous atmosphere, K(R ) is the dimensionless one-way irradiance extinction of the atmosphere at range 
R [m ], P(/? )(m - V - ‘] is the atmospheric aerosol backscatter coefficient, P L (t) is the pulse power profile 
of the laser, and C(R) is the coherent responsivity of the Doppler lidar flidar system dependence). The 
tenn P L (t-2R !c )[AT (R )] 2 B(/? )C (R ) describes the range weighting for the $NR. For many cases of 
interest, the terms [K(R )] 2 p(/? )C (R ) are constant over the range gate and the range weighting function is 
given by the term P L (t-2R !c ) The FWHM A r [m ] of the range weighting function for the Gaussian pulse 
is then 

Ar=cAt/2 (12) 

A simulation of Doppler Lidar data (Appendix A) and the range weighting function P L (t-2R /c ) is shown 
in Fig. 1 for typical parameters of a 2 micron lidar. At time 2psec the center of the pulse is located at 
range 300m with a width Ar . At a time 4 ps later the pulse is situated at a range of 900 m . The range 
resolution AR [m ] for an observation time MT$ is defined as 

AR=Ar+Ap (13) 

where Ap [m ] is the distance the range weighting function travels during the observation time, i.e., 

A p=MT s c/2. (14) 

The temporal scale of the modulation of the signal corresponds to the time required for the pulse to travel 
a distance Ar , after which time a new collection of independent atmospheric scatterers are illuminated by 
the pulse (see Fig. 1). 
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The periodogram of the data is defined as [Doviak and Zmic, 1984, 1993; Marple, 1987] 

T M — \ 

P(m y=-TT I £ z*exp(-2ju*m/A/)| 2 (15) 

M kM) 

and is shown in Fig. 1. This is an estimate for the spectrum at frequencies f =mAf where Af=l/(MT S ) is 
the frequency resolution. The maximum frequency that can be observed without aliasing is the Nyquist 
frequency F N =M A f =IIT S [for real data F N =l/(2T s )]. The frequency domain estimation of the mean fre- 
quency / involves extracting the location of the spectral feature of the random periodogram coefficients 
(or any spectral domain estimate) in the presence of random noise. It is reasonable to expect that the per- 
formance of a frequency domain estimator depends on the statistics of the signal and noise around the 
spectral feature. If one would increase the sampling interval T s by a factor of two and also reduce the 
noise bandwidth by a factor of two to match the new Nyquist frequency, the performance of the spectral 
domain mean frequency estimator would be unchanged because the statistics of the periodogram around 
the spectral feature is unchanged. However, the frequency range that can be searched has been reduced by 
a factor of two and the SNR has increased by a factor of two. Therefore, the SNR is not a convenient 
parameter to characterize the signal strength. 

A more useful parameter to describe the signal power for mean frequency estimation is the ratio of 
the total signal energy per estimate E =SMT s to the spectral noise level Nq=NT s , i.e., 

O=E/N 0 =SNR M (16) 

(Note that for real data, N q=2NT^ and <S>=SNR M/2). This is the basic parameter used in communication 
and radar for matched filter estimation of known signals in uncorrelated noise [Helstrom, 1968; Van 
Trees, 1968]. For Doppler lidar with quantum limited detection (the additive noise is dominated by the 
local oscillator shot noise), this parameter has a convenient physical interpretation. It is common practice 
to set the sampling rate MT S equal to the noise bandwidth B . The SNR Eq. (10) can then be written as 

SNR(ty=P R (tyr\ H (t)T s (17) 

where P R (t }=Po (z )/(h v)[photo -electrons I sec ] is the detection rate of photo-electrons. Therefore 

&*\hPr M T s . (18) 

Since P R M T s is the number of photo-electrons generated per observation time MT S , O is the effective 
number of photo-electrons coherently detected per observation time. For Doppler radar, this interpretation 
is not valid because quantum limited detection is not achieved. However, the interpretation of Eq. (16) is 
meaningful for Doppler radar. 

The spectral width w is commonly normalized by the Nyquist frequency F N =l/T s . For fixed range 
resolution (fixed A/7 s ), this parameter changes with the sampling interval. A more convenient parameter 
is the ratio of the spectral width w to the frequency resolution Af , 



0=w/A/ =wMT s 

(19) 

This parameter can be written as 

rt k Ar 

(20) 

where 

A/Tc 



'|s 

II 

£ 

(21) 


is defined as the number of independent samples of the signal in the observation time MT S . The parame- 
ter Q is also proportional to Ap/Ar, the number of pulse widths Ar contained in the distance Ap that the 
pulse travels during the observation time, if there are no velocity fluctuations over the sensing volume of 
the pulse that would increase the spectral width w . Using the parameters of Figure 1 we have M t =4.44, 
which is approximately the number of modulated events displayed in the lidar signal. For Doppler Lidar 
observations of uniform wind fields using a Gaussian pulse, the choice of one independent sample per 
observation interval (A//=l) produces 0=0.450 and Ap/Ar =2. Zmic (1979) investigated the number of 
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independent samples per estimate in terms of the reduction of the variance of the estimates for signal 
power. Zmic’s effective number of independent samples is also proportional to £2. The definition of Eq. 
(21) is chosen as a visual description of the data and as a simple connection to the covariance function. 
Since the choice of a definition for the number of independent samples per estimate is arbitrary, the 
parameter £2 will be used to represent this basic physical quantity. 

The covariance of the data [Eqs. (4), (5)] can be written in terms of the new variables <I> and £2, i.e., 

R k = ^-exp[-2Tt z Q 2 k 2 /M 2 +2Ki—£lk ]+5* . (22) 

M w 

The best performance of any unbiased estimator of the signal parameters is given by the Cramer-Rao 
Bound (CRB). For fixed O and £2, the CRB ay (Frehlich, 1993a) for estimation of the mean frequency 
normalized by the spectral width w as a function of M are shown in Fig. 2. [The normalized variable 
ay Iw is suggested by the functional form of Eq. (22)]. For sufficiently large M , the CRB for the normal- 
ized mean frequency is independent of M for all the parameters £2. For larger £2, the CRB becomes 
independent of M when the extent of the signal spectrum is less than the Nyquist frequency, i.e., when 
wT s <\/6 or A/>6Q. This condition is satisfied for any realistic problem. Since Of/w is independent of 
M , the CRB for mean frequency Of is independent of M for fixed w. This implies that the total observa- 
tion time per estimate MT S and range resolution is also fixed because Ci=wMT s is fixed. The CRB for 
mean frequency Of is independent of M for fixed observation time MT S ( fixed range resolution) and 
fixed w. This result was also obtained by Rye and Hardesty, (1993a) using approximate empirical 
models for Levin’s approximate CRB. The exact CRB verifies the accuracy of this statement. The CRB 
for estimation of the average signal power S and the CRB for the estimation of the spectral width w are 
also shown in Fig. 2. For both cases, the CRB is independent of M, provided that the spectral width is 
less than the Nyquist frequency. Similar results are produced for other values of <D. This result is excel- 
lent motivation for the choice of O and £2 as the basic system parameters. Accurate calculations of the 
CRB’s of the signal parameters for very large M can be made using the more efficient calculation at 
smaller Af. For fixed MT S (fixed range resolution), the frequency resolution A/ for the periodogram is 
also fixed. For fixed <D, £2, and MT S , the satisfies of the periodogram around the spectral peak will be 
unchanged because the shape of the signal peak compared to the frequency resolution Af is unchanged. 
This provides an intuitive explanation for the invariance of the CRB’s for the estimation of the normal- 
ized signal parameters. The low SNR analytic approximation to the CRB for mean frequency can be 
written as Of/w 2 = 4fnW<P 1 [Zmic, 1979, Eq. (A29)], which is independent of M. However, the high 
SNR analytic approximation to the CRB becomes oflw 2 = 12£2 2 /M 3 [Zmic, 1979, Eq. (A28)], which is 
not independent of M . It has been shown that the high SNR approximation to the CRB is not valid 
[Frehlich, 1993a]. 

The radial component v[m/s] of the velocity of the atmospheric scatterers is related to the mean- 
frequency / by \=Xf /2 and the velocity estimation error o v [m Is ]=XOfl2, where X[m ] is the operating 
wavelength. The velocity error is proportional to X if oy is independent of X, i.e., if <I> and £2 are indepen- 
dent of X and MT S (range resolution) is fixed. For Doppler Lidar O is independent of X for fixed laser 
power P L if the backscatter coefficient |3«X -1 «v [Eqs. (11), (16)]. £2 is independent of X if w is indepen- 
dent of X, i.e., there are negligible velocity fluctuations over the sensing volume of the pulse and w is 
determined by the transmitted pulse shape [Frehlich, 1993b]. Velocity error is then proportional to 
wavelength for fixed MT$ (fixed range resolution), fixed w, fixed power transmitted, and This is 

a convenient benchmark for system design and for some atmospheric conditions [Srivastava et al. 

1992]. Rye and Hardesty (1993a) found that Levin’s approximate CRB for mean frequency estimation is 
proportional to X m for fixed range resolution, fixed signal power, fixed maximum velocity, and fixed wT$ 
(the spectral width is determined by the velocity fluctuations over the sensing volume of the pulse). 

3. Mean Frequency Estimators 

For Doppler lidar with short pulses, the spectral width of the signal w will be determined by the 
pulse and therefore known a priori [Frehlich, 1993b]. If a long pulse is transmitted, the effects of wind 
shear and wind turbulence produce non-stationary data [Frehlich, 1993b]. For non-stationary data, spec- 
tral domain and covariance based estimators are not rigorously defined and the Maximum Likelihood 
estimators will produce the best performance. For a typical experiment, many pulses will be transmitted 
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and the SNR will also be known if the atmosphere is stationary over the total observation time. The 
mean- frequency estimators considered here will assume knowledge of the SNR and the spectral width w . 
For typical Doppler radar operation, the spectral width is determined by the wind fluctuations in the sens- 
ing volume, and must also be estimated. The pulse-pair algorithm does not require a priori knowledge. 
All the other algorithms discussed in this section require knowledge of the SNR and spectral width w . 


Pulse-Pair (PP) 

A numerically efficient algorithm for mean frequency estimation is the Pulse-Pair (PP) estimator 
given by [Miller and Rochwarger, 1972; Zmic, 1979] 

/-“=-««[*,] (23) 

2717,5 


where 



1 

M-k 


M-\ -k 
£ z,+*z* 


1=0 


is an unbiased estimate for R k , the covariance at lag kT s with a fixed length of data. 


(24) 


Maximum Likelihood Estimator (ML) 

The Maximum Likelihood (ML) estimator is the value of the unknown parameters that maximizes 
the log-likelihood function [Helstrom, 1968, Van Trees, 1968], given by 

L(f ,S ,w) = =-z* r R- 1 (/ r ^.H')z-ln[|R(A^,w)|]-A/ln(7t) (25) 

where | R | denotes the determinant of the matrix R. For the signal model considered here | R (f £ ,w ) | is 
independent of / . The simplest ML estimator for mean frequency / assumes knowledge of the signal 
power S and the spectral width w . For the covariance model of Eq. (4) [Zmic, 1979], 

Rj l (fJS,w) = D u (S ,w )cxp[2iufT s (k -/)] (26) 

where D(S ,w )=R _1 (/ =0,5 ,w ). In this case, the ML estimator for / is the value / that maximizes 

L i= £ d m cos(2nT s mf ) (27) 

m=0 


where 


M-m - 1 

d m = £ ZkZk+ n D kMn (S,w) (28) 

*=0 

The ML estimator can be efficiently calculated using the Fast Fourier Transform (FFT) algorithm of 
length M s >-M at discrete frequencies / =kfT s /M s . Note that the ML estimator is not given in terms of 
an estimate for the covariance. 


Periodogram Maximum Likelihood Estimator (PML) 

In the limit of large M , the DFT and periodogram coefficients P(m) become mutually uncorrelated 
and the Maximum Likelihood estimator for the mean frequency is the value of / that maximizes the log- 
likelihood function 


M- 


P(m) 


L(/^,w) = -£ 


(29) 


where [Doviak and Zmic, 1984, 1993, Chap. 5] 

P (m) = <P (m)> = T s X Pi exp(-2ruT s mllM) . (30) 

/=-(«- 1) 

[Levin (1965) and Rye and Hardesty (1993a, b) used this algorithm with the spectrum substituted for 
P (m J ,S ,w ). For small fl and large d>, this approximation results in poorer performance, e.g., for 0=0.5 
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and <t>>100 the velocity error is 15% larger than the PML estimator]. The PML estimator is numerically 
more efficient when computed in the time domain. If the DFT coefficients are mutually uncorrelated and 
statistically independent, then the inverse of the covariance matrix R is Toeplitz and has the form of Eq. 
(26) where 


D u (S,w) 


1 M ~ x exp[2 Kim(Jc -Q/M] 
M P(mJ =0,5 ,w ) 


(31) 


Substituting this result in Eq. (28) produces the periodogram ML estimator (PML), which is also 
equivalent to implementing Eq. (29) using the DFT algorithm. This estimator can be written in terms of 
estimates for the covariance R k because the matrix R " 1 is Toeplitz [See Eq. (28)]. 


Time Series Model Estimators (MVYW) 

Many spectral domain estimators have been proposed using models for the time series. These algo- 
rithms are discussed by Marple (1987) and computer programs provided. The regressive spectral estima- 
tor is given by 

7 — 02 ) 

1 1+S>* exp(-2r lifkTs ) | 2 

*=i 

where the coefficients a k and p*, are determined from the estimates of the covariance of the data. The 
order p is an input parameter that must be selected. 

Another time series algorithm introduced by Capon (1969) has been advocated by Anderson (1991) 
for Doppler lidar parameter estimation. Marple (1987) calls this algorithm a Minimum Variance (MV) 
estimator which can be written as 

1 p 1 

= V (33) 

Puv(pJ) &PAR(kf) 

For the spectral domain estimators, the mean-frequency / is determined by the maximum value of the 
spectral estimator. We have investigated the performance of the AR spectral estimator and the MV spec- 
tral estimator using the Yule-Walker and Burg method for estimating the regressive coefficients and using 
the order p that produced the best performance [smallest standard deviation of the good estimates (see 
Section 4)] for fixed (<h,Q M ). For large <D, the best order for the AR and MV estimators was p =1 , which 
corresponds to the PP estimator [Rye and Hardesty, 1993b]. The MV estimator with the Yule-Walker 
method (MVYW) produced the best overall performance. 

4. Performance of Mean Frequency Estimators 

The most useful description of the performance of an estimator is it’s Probability Density Function 
(PDF). For fixed T s and in the limit of large M , if the ML estimator is unbiased, it approaches the ideal 
performance of the CRB and t’s PDF approaches a Gaussian distribution. The periodogram (or any spec- 
tral estimator) of the random signal is random. For small d> or SNR , there will be realizations where the 
spectral feature of the signal is buried in the noise. Then, the spectral domain estimators for mean fre- 
quency will choose the largest noise feature that mimics the signal. This behavior was first recognized 
and analyzed by Shirakawa and Zmic, (1983) and Hardesty, (1986). This qualitative description of a 
clump of localized good estimates around the true mean frequency sitting on a pedestal of uniformly dis- 
tributed "bad" estimates has been observed in many other cases [Mahapatra and Zmic, 1983; May and 
Strauch, 1989; Anderson, 1991; Rye and Hardesty, 1993a], The PDF of the mean frequency estimators/ 
is therefore modeled as a Gaussian PDF centered on the true mean frequency / and a uniformly distri- 
buted component of bad estimates over the frequency range (0/^=1 /Ts), which is written as 

PDF (f)=bT s + —. ^- expH/-/ ) 2 /(2g 2 )] (34) 

yflKg 

where b is the fraction of bad estimates and g is the standard deviation of the good estimates (see also 
Rye and Hardesty, 1993c). The PDF of the mean- frequency estimators is estimated by it’s histogram with 
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bin size A/. The parameters of the model PDF (b ,g) are determined by minimizing the mean-square- 
difference between the histogram and the predicted histogram based on the model PDF. An example of 
the histograms for 10000 estimates of the ML estimator and the PP estimator, the predicted histogram of 
the best fit model parameters, and the model PDF are shown in Fig. 3. The ML estimator produces good 
agreement with the model PDF. The PP estimator is a poorer fit to the model PDF when <t><100. This is 
due to the nonlinear interaction of the signal and noise. 


The comparison of estimator performance with the CRB is complicated by the presence of the uni- 
form distribution of bad estimates. If the model PDF is a good approximation then it provides a complete 
statistical description of an estimator. For example, the ensemble average of the mean frequency estima- 
tor/ is 


l/Ti 


</> = j f PDF (J )df = 
o 


b 

27V 




(35) 


which indicates a bias from the true mean-frequency /. The mean-square-error of the estimator is 
defined as 

mS h 

<(f~f) 2 > = °e 2 = j (/-f) 2 PDF(f)df = ~ 1 (l-3fT s +3f 2 T s 2 }Hl-b)g 2 (36) 

o 37V 


which has a minimum for /7V =0.5. The CRB for a biased estimator is given by [Helstrom, 1968] 

a 

«/-/) 2 > S = (l-b) 2 o f 2 (37) 

where ay 2 is the CRB for an unbiased estimate of / . The CRB for a biased estimate depends on the frac- 
tion of bad estimates b which in turn depends on the estimator/. A comparison of different estimators 
that includes the bias can be related to the CRB for unbiased estimates by dividing Eqs. (36) and (37) by 
(1 -b) 2 . From the modified Eq. (36), the adjusted standard deviation of the estimator error becomes 


o c =(V(l -b) (38) 

which is to be compared with the CRB Of for unbiased estimators. 

An empirical model for the fraction of bad estimates b as a function of <I> for fixed fl and M is 

6(<t)=[l+(0/6o) a r. (39) 

In the previous section, it was shown that the CRB for estimation of the mean frequency / normalized by 
the spectral width w was independent of the number of observation points M for fixed O and Cl. This 
suggests the use of the parameter gs=g/w • An empirical functional form for this parameter is 

gN (<*>) = Xt l-KWg (40) 

The simulation results shown in the following Figures were generated by averaging 10 simulations of 
2000 realizations each for large <E>. As the number of good estimates declined at small <D, the number of 
realizations were increased to maintain the number of good estimates at approximately 2000. The stan- 
dard deviation of the 10 estimates for b and gs-g /w were used to obtain the chi-squared best fit [Press et 
al., 1986] for the parameters of the empirical models Eqs. (39) and (40). 

The performance of the Maximum Likelihood estimator is shown in Fig. 4 for Q=0.5. The simula- 
tion results agree well with the best fit models Eqs. (39) and (40). For large <I>, the bias of the estimator is 
small and the standard deviation of the good estimates g approaches the CRB. For small <£, the com- 
parison of the standard deviation c c of the biased estimate to the CRB that includes the effect of the bias 
is poor. This is a typical result for all the estimators we investigated. The CRB with the bias included [Eq. 
(37)] is a poor indicator of estimator performance. The standard deviation of the good estimates g and 
the fraction of the bad estimates b provide a complete description of the performance of mean-frequency 
estimators when the model PDF defined by Eq. (34) is valid. For small O, the standard deviation of the 
good estimates g is less than the CRB. This does not contradict the meaning of the CRB because in this 
regime the estimates are biased and are not bounded by the CRB for unbiased estimates. 
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The performance of the periodogram maximum likelihood (PML) estimator is shown in Fig. 5 for 
£2=0.5. The simulation results agree well with the best fit models Eqs. (39) and (40). When the estimates 
are unbiased (<t>>~100), the performance is bounded by the periodogram CRB. 

The performance of all the mean-frequency estimators considered here are shown in Fig. 6 for 
£2=0.5 (typical Doppler lidar case) and in Fig. 7 for 0=2.0 (typical Doppler radar case). The parameters 
of the empirical models Eqs (39) and (40) are given in Table 1. All the spectral domain and covariance 
based estimators are bounded by the periodogram CRB for <t>>~100 (unbiased estimates). For the best 
estimators (ML,PML,AR,MV), the standard deviation of the good estimates g approaches w for small O 
and approaches a fraction of w (=1/10) for large $>. This general feature produces a useful description of 
mean- frequency accuracy in terms of the basic parameters $ and Q. The dependence of b with <t> is 
similar for the best estimators. The PDF of the PP estimator is a poor fit to the empirical model Eq. (34) 
for <D<~100 (see Fig. 3) and results in deviation from the other estimators. For 0=0.5 and large O, more 
than a factor of two improvement in mean- frequency accuracy is realized with the Maximum-Likelihood 
estimator compared with the best spectral -domain and covariance based estimators. 

The performance of the Maximum Likelihood estimator for A/ =32,64,128 is shown in Fig. 8 for 
£ 2 = 0.5 and in Fig. 9 for £2=2. The standard deviation of the good estimates has little dependence on the 
total number of data points M for <I>>10, i.e., the mean- frequency error is independent of the maximum 
observable frequency F N for fixed range resolution and fixed signal energy per estimate. For fixed O, and 
increasing M , the fraction of bad estimates b increases because there are more noise spikes that may be 
larger than the true mean-frequency signal peak of the log-likelihood function. 


5. Tradeoff Between Pulse Energy and Number of Pulses 

Rye and Hardesty (1993b) posed an important question: is it better to use many pulses with low 
energy or one pulse with the same total energy? For constant laser power and for the case where accumu- 
lated data from many pulses approaches Levin’s CRB, they showed that many low energy pulses were 
better than one high energy pulse for narrow spectral width. This regime corresponds to large <t> where 
there are no bad estimates. The performance of mean-frequency estimators for K pulses per estimate can 
be compared to the single pulse case for fixed <D. The data for each pulse is assumed statistically 
independent. This is a good assumption for Doppler Lidar because the signal correlation time is of the 
order of lps and the interval between pulses is larger than 10msec . 

The Likelihood Function L K for K statistically independent data vectors z„ and fixed range resolu- 
tion (fixed MT S ) is [Scharf, 1991] 

K 

L k (f ,S K ,w ) = -K\n[\R(f JS K ,w)\] - KM\n(K}-'Z z « T ' S K z « (41) 

rt=l 

where z T denotes the transpose of z and S k =S/K is the signal power for each pulse. The CRB for this 
data set of K pulses is 

ofc=o f \S K )/K (42) 

where Or 2 (Sk ) is the CRB for one pulse with signal power S*. The CRB for the mean- frequency as a 
function of the total number of pulses but with fixed total energy 4> is shown in Fig. 10. Multiple pulses 
are better for <D>~10 and smaller fi. The performance of the ML estimator [the value of / that maxim- 
izes Eq. (41)] for different K as a function of O is shown in Fig. 1 1 for 0=0.5 and in Fig. 12 for 0=2. For 
small <t>, the best performance is obtained by using one pulse. For large O, performance is improved by 
dividing the total laser energy among many pulses. The effects are most pronounced for small O. For 
small <I>, concentrating all the energy into one pulse produces the optimal signal for identifying the good 
estimates. For large O, there is no difficult identifying the good estimates and more pulses result in a 
decrease in the mean-frequency error because each pulse produces a statistically independent estimate and 
the error scales as K- V2 o f (S K ). Since o f (S K ) changes slowly compared to AT 1 ' 2 , performance improves 
with more pulses. The CRB (see Fig. 10) provides an effective prediction of the behavior of the standard 
deviation of the good estimates g for all <1>. 

6. Summary 
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For the signal model of Eqs. (2), (4), and (5), the perfoimance of mean-frequency estimators is con- 
veniently presented in terms of two basic parameters: O the ratio of the signal energy per estimate to the 
spectral noise level and Q which is proportional to the number of independent samples per estimate. For 
Doppler lidar, <t> is the number of coherently detected photo-electrons per estimate. For fixed <t> and Q, 
the CRB’s for unbiased estimates of the signal parameters are independent of the number of samples per 
estimate M (Fig. 2). The PDF of a mean-frequency estimator is the best description of performance. For 
the better estimators, the PDF can be accurately modeled as a Gaussian distribution function around the 
true mean frequency and a uniformly distributed component of bad estimates (Fig. 3). This model 
requires two parameters: the standard deviation of the good estimates g and the fraction of the estimates 
b that are bad. The empirical functions Eqs. (39) and (40) provide a good approximation for g and b as a 
function of O for fixed Q and \ M . The parameters of these empirical models are listed in Table 1 for some 
common mean-frequency estimators. The fraction of bad estimates can be reduced by neglecting the esti- 
mates that have low signal power (Rye and Hardesty, 1993c). Adaptive mean-frequency estimators have 
also been proposed to reduce the detrimental effects of the fades in signal power [Dabas et al., 1993], 

For large <t>: the estimates are unbiased, the standard deviation of the spectral domain and covari- 
ance based estimates are bounded by the periodogram CRB (Figs. 6,7), and the standard deviation of the 
ML estimator approaches the exact CRB (Fig. 4). For smaller Q, the ML estimator has more than a factor 
of two improvement in performance compared to previous estimators (Fig. 7). For <D<~100, all the esti- 
mators are biased due to the uniformly distributed bad estimates. For fixed M , the behavior of the fraction 
of bad estimates b is almost a universal function of O with a weak dependence on Q (Figs. 6-9), i.e., the 
better estimators (ML,PML,AR,MV) perform equally well in identifying the cases where the random sig- 
nal is not buried in the random noise. For fixed 0> and fl, as M increases (Figs. 8,9), the fraction of bad 
estimates increases because the random spectral feature of the signal must be identified among more noise 
spikes. For <t>>100 and fixed Q, the standard deviation of the good estimates normalized by w is essen- 
tially independent of M and for fixed range resolution and fixed signal energy, the standard deviation of 
the mean frequency estimators are independent of the maximum observable frequency F N . For fixed 
range resolution, fixed spectral width, fixed transmitted pulse energy, and a backscatter coefficient p^X -1 , 
the velocity error is proportional to the operating wavelength X. The spectral width of the signal is con- 
stant with operating wavelength when the wind fluctuations are negligible over the sensing volume of the 
pulse. This is true for Doppler lidar with short pulse lengths and typical atmospheric conditions. For typi- 
cal Doppler radar conditions, the spectral width of the signal is determined by the wind fluctuations over 
the sensing volume of the pulse. Then, the velocity error is approximately proportional to X ,/2 for fixed 
range resolution and fixed signal energy (Rye and Hardesty, 1993a). 

For Doppler lidar and for large O, it is better to transmit many low energy pulses instead of one 
pulse with the same total energy, if the parameters of the process do not change during the measurement 
interval. For small O, the converse is true (Figs. 11,12). 

The behavior of the performance of mean-frequency estimators as a function of the two basic 
parameters O and simplify Doppler Lidar design. In many cases, the range resolution of the measure- 
ment is determined by the application. Assuming one independent velocity measurement per observation 
time and small fluctuations of the wind field over the sensing volume of the pulse (w is determined by the 
transmitted pulse profile and fW).45) fixes the pulse length and observation time per estimate MT S . The 
maximum search velocity v mix =X/(2T i j) determines T$ and M . The velocity accuracy Oy=Xoy/2 and the 
fraction of bad estimates b is then a function of one parameter <t>. To produce estimates with a small 
number of bad estimates, <&>100, which determines the required pulse energy for a specified backscatter 
coefficient. 


Acknowledgement. The authors acknowledges useful discussions with M. J. Kavaya, D. S. Zrnic, and B. 
I. Rye. This work was supported by the National Science Foundation and the National Aeronautics and 
Space Administration. 

APPENDIX A: Simulation of Doppler Radar/Lidar Data 


Revised September 22, 1993 


11 


Complex Doppler radar/lidar data with a specified covariance can be generated with a simple com- 
puter algorithm. The complex Gaussian data vector z is produced by 

Ms - 1 

z k = 2 y* exp(2iu'/i* /Ms ) (Al) 

n = 0 

where y n are zero mean complex Gaussian random variables with 

<yny* m > = Y n 5 n - m (A2) 


and 

i Ms - 1 

V„ = — Y R k <ixp(-2nink/M s ) (A3) 

*=o 

where /?* is the desired covariance of the stationary data. The random vector z then has the properties 


<Z*Z* > = /?*_/ = <Z k ZM s -l> 


(A4) 


and the first points of the simulated data are correlated with the end points. Therefore, a simulation of 
length M S =2M is required to generate an independent realization of length M . The length M s must be 
long enough to ensure that all Y n are positive [see Eq. (A3)]. Zmic (1975) and Sirmans, and Bumgarner, 
(1975) used this algorithm with Y n chosen to produce exact sampling of the spectrum instead of an exact 
sampling of the covariance. For Doppler lidar, the covariance of the data is specified by the transmitted 
pulse and lidar system parameters [Frehlich, 1993b]. 
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Figure Captions 

Fig.l Simulation of Doppler lidar data for a 2 pm lidar [Henderson et al, 1993], with /V=1.0, SNR= 10, 
f=5MHz, w=Q.2MHz , T s =20nsec . The range weighting function P L {t-2R!c) is shown for the 
signal at 2 ps and 6pj . The FWHM A r of the pulse sensing volume, the distance A p that the pulse 
moves in 4pj , the signal correlation time At, and the periodogram Eq. (15) are also shown. 

Fig.2 The normalized CRB ’s oy for estimation of mean frequency / , o s for estimating signal power S , 

and G w for estimating spectral width w as a function of M for 0=10 and Q = 0.2 ( ), 0.5 (....), 

1.0 (- - -), 2.0 (. - . -), and 5.0 . -). 

Fig.3 The histograms (+) from 10000 ML and PP estimates of mean frequency/ =k/Af where kf is the 
mean-frequency index, 0=10, 0=0.5, and A/ =32. The best fit empirical model Eq. (34) for the 
PDF is indicated by ( ) and the histogram determined from the best fit PDF is indicated (O). 

Fig.4 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad" estimates as a function of O for 0=0.5 and M =32. The results of the simulation are indi- 
cated by (O) and the best fit empirical models Eqs. (39) and (40) are indicated by ( ). The 

exact CRB is (....) and the empirical model corrected for bias by Eq. (38) is indicated by ( — ). 

Fig.5 The standard-deviation g of the "good" PML estimates for mean-frequency and the fraction b of 
"bad" estimates as a function of O for 0=0.5 and M= 32. The results of the simulation are indi- 
cated by (O) and the best fit empirical models Eqs. (39) and (40) are indicated by ( ). The 

exact CRB is (....) and the periodogram CRB is indicated by ( — ) [Frehlich, 1993a]. 

Fig.6 The standard-deviation g of the "good” estimates for mean-frequency and the fraction b of "bad" 
estimates as a function of <1> for 0=0.5 and Af=32. The results of the simulation are given by the 

best fit empirical models Eqs. (39) and (40) for the ML ( ), PML (....), MVYW, ( — ), and PP 

(. - • -)• 

Fig.7 The standard-deviation g of the "good" estimates for mean-frequency and the fraction b of "bad" 
estimates as a function of O for 0=2.0 and M =32. The results of the simulation are given by the 

best fit empirical models Eqs. (39) and (40) for the ML ( ), PML (....), MVYW, ( — ), and PP 

(.-.-). 

Fig.8 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad" estimates as a function of <I> for 0=0.5. The results of the simulation are given by the best fit 
empirical models Eqs. (39) and (40) for M = 32 ( ), 64 (....), and 128 ( — ). 

Fig.9 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad" estimates as a function of d> for 0=2.0. The results of the simulation are given by the best fit 
empirical models Eqs. (39) and (40) forM = 32 ( ), 64 (....), and 128 ( — ). 

Fig. 10 The CRB for estimation of mean frequency for multiple pulses with fixed total energy (fixed <I>) 
and 0= 0.2 ( ). 0.5 (....), 1.0 (- - -), 2.0 (. - . -), and 5.0 (.. - .. -). 

Fig. 11 The standard-deviation g of the "good" ML estimates for mean-frequency for Doppler lidar and 
the fraction b of "bad" estimates as a function of <I> for different number of pulses per estimate 
with the total signal energy per estimate fixed and 0=0.5. The results of the simulation are given 

by the best fit empirical models Eqs. (39) and (40) for 1 pulse ( ), 5 pulses (....), and 10 pulses 

(---). 


Revised September 22, 1993 



15 


Fig. 12 The standard-deviation g of the "good" ML estimates for mean-frequency for Doppler lidar and 
the fraction b of "bad" estimates as a function of O for different number of pulses per estimate 
with the total signal energy per estimate fixed and £2=2.0. The results of the simulation are given 

by the best fit empirical models Eqs. (39) and (40) for 1 pulse ( ), 5 pulses (....), and 10 pulses 

(- - -). 
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Lidar Signal Range Weighting 


Fig.l 


Simulation of Doppler lidar data for a 2 pm lidar [Henderson et al, 1993], with N =1.0, SNR =10, 
f =SMHz, w=C.2MHz, T s =2Dnsec . The range weighting function P L (t-2Rlc) is shown for the 
signal at 2\ls and 6|is . The FWHM Ar of the pulse sensing volume, the distance Ap that the pulse 
moves in , the signal correlation time At, and the periodogram Eq. (15) are also shown. 
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Figure 2 Frehlich JTech No. 93052002 


Fig.2 The normalized CRB’s Of for estimation of mean frequency / , o s for estimating signal power S , 

and o w for estimating spectral width w as a function of M for <t>=10 and Q = 0.2 ( ), 0.5 (....), 

1.0 (---), 2.0 and 5.0 
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Figure 3 Frehlich JTech No. 93052002 


Fig.3 The histograms (+) from 10000 ML and PP estimates of mean frequency/ =k/Af where kf is the 
mean-frequency index, <I>=10, ft=0.5, and A/ =32. The best fit empirical model Eq. (34) for the 
PDF is indicated by ( ) and the histogram determined from the best fit PDF is indicated (O). 
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Figure 4 Frehlich JTech No. 93052002 


Fig.4 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad” estimates as a function of <t> for £>=0.5 and Af=32. The results of die simulation are indi- 
cated by (O) and the best fit empirical models Eqs. (39) and (40) are indicated by ( ). The 

exact CRB is (....) and the empirical model corrected for bias by Eq. (38) is indicated by ( — ). 
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Figure 5 Frehlich JTech No. 93052002 


Fig.5 The standard-deviation g of the "good" PML estimates for mean- frequency and the fraction b of 
"bad" estimates as a function of <J> for 0=0.5 and M =32. The results of the simulation are indi- 
cated by (O) and the best fit empirical models Eqs. (39) and (40) are indicated by ( ). The 

exact CRB is (....) and the periodogram CRB is indicated by ( — ) [Frehlich, 1993a]. 





Figure 6 Frehlich JTech No. 93052002 


Fig.6 The standard-deviation g of the "good" estimates for mean- frequency and the fraction b of "bad" 
estimates as a function of <I> for £>=0.5 and M =32. The results of the simulation are given by the 
best fit empirical models Eqs. (39) and (40) for the ML ( ), PML (....), MVYW, ( — ), and PP 

(. - • -)• 
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Figure 7 Frehlich JTech No. 93052002 


Fig.7 The standard-deviation g of the "good” estimates for mean-frequency and the fraction b of "bad” 
estimates as a function of <I> for te=2.0 and M =32. The results of the simulation are given by the 
best fit empirical models Eqs. (39) and (40) for the ML ( ), PML (....), MVYW, ( — ), and PP 

(•-.-)• 
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Figure 8 Frehlich JTech No. 93052002 


Fig.8 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad" estimates as a function of O for £2=0.5. The results of the simulation arc given by the best fit 
empirical models Eqs. (39) and (40) for M = 32 ( ), 64 (....), and 128 ( — ). 
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Figure 9 Frehlich JTech No. 93052002 


Fig.9 The standard-deviation g of the "good" ML estimates for mean-frequency and the fraction b of 
"bad” estimates as a function of <I> for £1=2.0. The results of the simulation are given by the best fit 
empirical models Eqs. (39) and (40) for M = 32 ( ), 64 (....), and 128 ( — ). 




Figure 10 Frehlich JTech No. 93052002 


Fig. 10 The CRB for estimation of mean frequency for multiple pulses with fixed total energy (fixed d>) 
and Cl= 0.2 ( ), 0.5 (....), 1.0 (- - -), 2.0 (. - . -), and 5.0 (.. - .. -). 
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Figure 11 Frehlich JTech No. 93052002 


Fig. 11 The standard-deviation g of the "good" ML estimates for mean-frequency for Doppler lidar and 
the fraction b of "bad" estimates as a function of <I> for different number of pulses per estimate 
with the total signal energy per estimate fixed and 0=0.5. The results of the simulation are given 
by the best fit empirical models Eqs. (39) and (40) for 1 pulse ( ), 5 pulses (....), and 10 pulses 
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Figure 12 Frehlich JTech No. 93052002 


Fig. 12 The standard-deviation g of the "good” ML estimates for mean-frequency for Doppler lidar and 
the fraction b of "bad” estimates as a function of O for different number of pulses per estimate 
with the total signal energy per estimate fixed and £1=2.0. The results of the simulation are given 

by the best fit empirical models Eqs. (39) and (40) for 1 pulse ( ), 5 pulses (....), and 10 pulses 

(- - -)• 
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ABSTRACT Calculations of the exact Cramer-Rao Bound (CRB) for unbiased estimates of the mean fre- 
quency, signal power, and spectral width of Doppler radar/lidar signals (a Gaussian random process) are 
presented. Approximate CRB’s are derived using the Discrete Fourier Transform (DFT). These approxi- 
mate results are equal to the exact CRB when the DFT coefficients are mutually uncorrelated. Previous 
high SNR limits for CRB’s are shown to be inaccurate because the discrete summations cannot be 
approximated with integration. The performance of an approximate maximum likelihood estimator for 
mean frequency approaches the exact CRB for moderate signal to noise ratio and moderate spectral 
width. 

EDICS - 2.1, 4.1, 4.9.1, 5.1, 5.2, 6.3.2 


I. INTRODUCTION 

The estimation of the mean frequency, signal power, or spectral width of radar, sonar, lidar, and 
Synthetic- Aperture-Radar (SAR) returns [1-28] is important for remote sensing. Many estimation algo- 
rithms have been investigated, including pulse pair [1 1,24,26], peak firequency[4], spectral domain [1,3- 
13,17,18], maximum entropy [8,9,17,18], and maximum likelihood [1,3,7,17,18]. The use of accumu- 
lated periodograms or correlograms for mean frequency estimation has been investigated by Rye and Har- 
desty [17,18]. Many spectral domain algorithms have been reviewed by Kay and Marple [7], 

Doppler radar data from atmospheric returns is sampled from the output of a complex receiver 
[29,30]. For Doppler radar the scattering targets are atmospheric refractive index irregularities, rain, 
snow, or insects. For Doppler Lidar the scattering targets are atmospheric aerosols. In both cases, the 
received signal is well represented as a Gaussian random process because the total backscattered field is 
the superposition of many randomly phased scattered fields. For Doppler radar, the time series is gen- 
erated by transmitting a series of pulses. For Doppler lidar, the time series is the return signal from a sin- 
gle transmitted pulse. The spectral width of a Doppler lidar signal is usually smaller than the spectral 
width of a Doppler radar signal for the same range resolution because the effects of velocity fluctuations 
over the sensing volume of the lidar pulse is smaller than for the radar case. The mean frequency of the 
return signal provides an estimate of the radial component of the velocity field averaged over the sensing 
volume of the pulse [2], 

The ideal performance of an estimator of the signal parameters is given by the Cramer-Rao Bound 
(CRB) [16,30,31]. The exact CRB has been determined for the cases of a sine wave with random ampli- 
tude and phase in white noise [26], and for an autoregressive signal model [16]. The exact CRB for 
parameter estimation of a Gaussian random process have been derived but not applied to the Doppler 
radar/lidar problem. A new, more efficient expression for the exact CRB for estimates of the mean fre- 
quency of a complex Gaussian signal added to complex Gaussian uncorrelated noise is derived. For this 
data model, it is shown that the exact CRB’s of all the signal parameters are independent of the mean fre- 
quency. Approximate calculations of CRB have been discussed by Whittle [21], Zmic [26], and Levin 
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[10]. These approximations arc essentially large sample approximations. Approximate CRB’s are 
developed by demonstrating that the CRB of a complex Gaussian process is unchanged by the Discrete 
Fourier Transform (DFT), a unitary transformation. Approximate CRB’s for the signal parameters are 
derived using the periodogram (squared modulus of the DFT) of the data. These periodogram CRB’s are 
equal to the exact CRB when the periodogram elements are mutually uncorrelated. The maximum corre- 
lation coefficient can be calculated thus providing an efficient test for the validity of the approximate 
CRB. The periodogram CRB’s and Levin’s approximate CRB’s are compared with the exact CRB’s for a 
complex Gaussian signal in complex Gaussian uncorrelated noise with a Gaussian covariance function for 
the signal. The performance of an approximate maximum likelihood estimator (ML) for the mean fre- 
quency / is compared to the CRB. The validity of previous analytic approximations for the CRB for 
mean frequency are determined. A summary concludes the paper. 


H. THEORETICAL DEVELOPMENT 


A. Exact Cramer-Rao Bounds 

Doppler radar signals are well represented as a zero-mean complex Gaussian random process 
[29,30] z* where z is a column vector of M samples, 

E[z k z,] = 0 (1) 


and £[.] denotes ensemble average. The statistics of z are completely determined by the covariance 
matrix [31,33] 

R(9) = £[z.z* r ] R*,(0) = £[z*z,*] (2) 

where 8 is the vector of parameters of the process, the row vector z r is the transpose of z, * denotes com- 
plex conjugation, and Ry denotes the elements of the covariance matrix R. 

The variance of an unbiased estimate ($, of the parameter 0, has a lower bound called the Cramer- 
Rao bound (CRB) which is given by [18,31] 

VAR [$,-] =E[(d 1 -e, ) 2 ] >/“(0) (3) 


where J“ is the i th element of the diagonal of the inverse of the Fisher information matrix J, with ele- 
ments given by [18,31] 


d 2 ln[ | R(8) | ] 
90 * 39 / 


+ £ 


,*T 


9 2 R- ] (8) 

90*90, 


z 


9 2 ln[ | R(8) | ] , , [- 3 2 R- 1 (0) ' 
90* 90, 30, 90, 


(4) 


where | R | denotes the determinant of the matrix R and rr< 


R 


denotes the trace of matrix R. This result 


is difficult to evaluate because the inverse of the matrix R is required before the derivatives are computed. 
However, for a zero mean complex Gaussian process [17,15] 


/*, = iHR ‘ 1 ( 0)-® D - 1 


r - V8 ,M 

30* - (8) 30, 


(5) 


Exact calculations of the CRB simplify for the common model of Doppler radar or Doppler lidar 
atmospheric signals given by 

z* = s k txp(2nifkT s ) + n k (6) 

where s k and n k are complex Gaussian random processes, / [Hz ] is the mean frequency of the signal, 
T s [ sec ] is the sampling interval, the noise samples «* are statistically independent of the signal s k , 

E[n k n t } = 0 £[**«,•] =N8 W , (7) 

N is the average noise power per sample and 8* is the Kronecker delta symbol. The statistics of s k are 
assumed to be independent of the mean frequency / . The covariance matrix R then becomes [26] 
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R(0) = DQ(4>) D* 


( 8 ) 


or 


*«(B) = Q U mxp[2KifT s (k-l)] 
where D is a diagonal matrix with elements 

D h t =exp(2jt ifkT s ) , 


(9) 

( 10 ) 


Q u are the elements of the covariance matrix Q, the covariance of s k +n k , and the system parameters are 
e =(M). The matrix D is a unitary matrix [DD* = I where I is the identity matrix]. The determinant of 
R is equal to the determinant of Q and is independent of the mean frequency / . 


The Fisher information matrix for elements involving the parameters + are given by 

jgai ] 


y a = ,Jd r'M 


srVgfffw * i 


( 11 ) 


since 


jDQDr'h^ 


( 12 ) 


for any matrices D and Q. These Fisher information matrix elements are independent of the mean fre- 
quency / . The same argument proves that the Fisher information matrix elements involving the mean fre- 
quency are also independent of the mean frequency. Therefore, the Fisher information matrix and the 
CRB’s are independent of the mean frequency/ for the signal model of Eqs. (6) and (7). 

A numerically more efficient expression for the exact CRB of estimators for the mean frequency is 
produced by writing the inverse of the covariance matrix as 

R- 1 (9) = D^- 1 (4>)D- 1 (13) 


or 

(R-')tf (0) = (£~V*)exp[2ru/r 5 (*-/)] (14) 

where (R -I ) w denotes the kl element of the matrix R" 1 . The Fisher information matrix element J u for 
the mean frequency/ becomes [substitute Eqs. (9) and (14) into Eq. (4)] 

j u = - 4 J c 2 r s 2 5: S(^-*) 2 G*/(S _1 )tt O 5 ) 

*=0 1=0 

This result can also be obtained from Eq. (5) by expanding the matrix operations and using the matrix 
identity QQ ^I. Eq. (5) has approximately 3M 3 +2Af 2 +Af floating point operations while Eq. (15) has 
approximately 2M 2 +M floating point operations. If the mean frequency / is the only unknown parame- 
ter or if the off diagonal elements of J are zero, the CRB for estimation of the mean frequency / is given 
by 

cj>\!J n . (16) 


B. Approximate Cramer-Rao Bounds 

The calculation of the CRB is numerically intensive, especially for large M . Asymptotic expres- 
sions valid for large M have been derived by Whittle [21], Levin [10], and Zmic [26]. The region of vali- 
dity of these results is difficult to determine [10,15] and depends on the parameters of the process as well 
as M. A new approximate CRB and a criteria for testing the validity of the approximation will be 
presented here, based on the past work of Levin and Whittle. The approximate CRB’s of Whittle, Levin, 
and Zmic can be derived from this approximate CRB, however, further approximations are required. 

The Discrete Fourier Transform (DFT) of the data vector z is given by 

u = U z (17) 
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where the elements of the unitary matrix U are 

U u =M- ll2 expi-2KikllM ) 


«* =M 1/2 £ zitx^-lnikllM) . (19) 

1=0 

Because this is a linear transformation, the vector u has a joint complex Gaussian probability density 
function which is a function of the covariance matrix R' 

R' = E [u u* T ] = E [U z z* T y 1 ] = U R U" 1 . (20) 

The determinant of the matrix R' is equal to the determinant of the matrix R. The Fisher information 
matrix for the parameters 0 in terms of the statistics of the new vector u is denoted as J' and given by [see 
Eq. (4) and (12)] 


J'h = 


5 2 ln[|R(0)|] 

50 * 50 / 




50*50/ 


The Fisher information matrix using the random vector u is equal to the Fisher information matrix using 
z. Therefore, the CRB’s are identical. 

The vector of DFT coefficients u contain the same information as the original time series z. 
Assume the random variables «* are mutually uncorrelated and statistically independent [10]. Then R' is 
a diagonal matrix with elements 

/?'*/= £[«*«/*] = />(*, 0 ) 8 W . ( 22 ) 

The function P ( k ,0) is the ensemble average of the periodogram 

P(k) = u k u k (23) 

and can be written as [28] 

P(*,0)= M £ /?„(©)( 1 - ~^)exp(-2ju*n/A/) (24) 

n — ( M - l ) M 

and /?*/ (©)=/?*_/ (0). The inverse covariance matrix R' -1 is a diagonal matrix with elements 

(R'- 1 )*/=P- I (i,0)8*_/ . (25) 

Then, the Fisher information matrix is given by [substitute Eqs. (22) and (25) into Eq. (5)] 

M-id(m,0*)rf(m,0/) 

J ki-Jki= L 2(m Q, < 26 > 

m=o P \mM) 


where 


5(m,0*) = 


5 P ( m , 0 ) 


For finite M , the accuracy of this approximation for the CRB is determined by the correlation coefficient 
p(* ,/) of the periodogram elements which can be written as [28] 


P( * ' l) P(k ,B)P (l ,0W 2 1 (G) 


sin[ft|w | (A: - l)!M] 
sin[7t(/ -k)/M] 


exp[-2 ni (k + 1 )n/M] | 2 . (28) 


The results of Whittle and Levin can be obtained from Eq. (26) by assuming the average periodo- 
gram can be approximated by the continuous spectrum and appro ximating the summation as an integral. 
Zmic [26] produced an approximate time-domain CRB for the estimation of the mean frequency / and 
the signal model Eqs. (6) and (7), by assuming R -1 , the inverse of the covariance matrix R, is asymptoti- 
cally Toeplitz (see ref. [26], Eqs. (A.21), (A.26)). This result can be derived from the approximate 
periodogram CRB [34], 
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UI. CALCULATIONS 


A. Numerical Calculations 

The CRB is calculated for the signal model of Eqs. (6) and (7), assuming stationary of the data 
(/?m(0)=/?*_^(0)] and the Gaussian covariance model of Zmic [26], i.e., 

R k =S exp[-27t 2 (w'&7s ) 2 +2nifkT s ]+N 6* (29) 


or 


Q k =S expl-h^iwkTs ) 2 ]+N 5* 


(30) 


where 5 is the signal power, / is the mean frequency or Doppler frequency, w the spectral width, N is 
the noise power, and SIN is the Signal-to-Noise Ratio (SNR). The system parameters are 
( 01 , 02 . 03 ) = (f ,S,w ). Then 

P(m,0) = N+S X (1 — ~~)exp[— 2jt 2 w 2 T5 2 n 2 ]exp[— (m . (31) 

hm-1) m M 

The calculations of the approximate CRB simplify if we write the mean frequency / as 


/ = {m F + e)A f 


(32) 


where 


A / = 


1 

MT S 


(33) 


is the frequency resolution of the periodogram, m F is the mean frequency index, and —1/2 <e< 1/2. Then 
P (m ,0) = N + S X O “ -m F - e)] . (34) 

n=-(W-l) M M 

The integer m F is the index of the maximum of the average periodogram and e is the offset in units of A f 
of the true mean frequency from the frequency m F Af . The elements of the Fisher information matrix for 
the approximate CRB can be written as a DFT. The function d(m ,© 2 ) is the coefficient of S in Eq. (34). 
The other two functions can be written as [see Eq. (31)] 




dP Qn ,0) 

d/ 


= 2niST s 


M-\ 


Z (1- 

n =-(M —1) 


M 

M 


)n 


expl-li^w 2 T s 2 n 2 ]exp[— ^^-(m -m F -e)] 

M 


(35) 


and 


e 3)= dP(k,Q) = ^ v 2 S wT s 2 x' (1- 


d(m 




(M-l) 


l«l 

M 


)n‘ 


2nin 


txp[-2T^w 2 T s 2 n 2 ] l txp[ — —(m -m F -e)] . 


M 


(36) 


The approximate CRB is a function of the frequency offset parameter e. For spectral width 
wMT s >0.5, there is little dependence on e [34], The mean frequency / is an unknown parameter which 
is usually random. Then e will be uniformly distributed on the interval (-1/2.1/2). To remove the depen- 
dence of e, the approximate CRB’s are calculated using Eqs. (4), (16), (26), (34), (35), (36) and the result 
averaged over the uniform distribution of e. We will call this approximate CRB the periodogram CRB. 

The correlation coefficient p (k ,/ ) is also a function of e. The average over e of the maximum corre- 
lation coefficient for all k and l is denoted by p max . If p max is small, the periodogram elements are 
approximately mutually uncorrelated and the periodogram ORB’s are a good approximation for the exact 
CRB’s of unbiased estimates of the parameters of the complex Gaussian random processes u and z. For 
all periodogram CRB calculations, the index m F is set to zero since the summation over index m in Eq. 
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(26) is independent of trip. The numerical calculations of the periodogram CRB are performed using the 
Fast Fourier Transform (FFT) algorithm for the DFT. For the exact CRB, the inverse of the Toeplitz 
matrix Q is calculated using Trench’s technique with Durbin’s algorithm [35]. 


Zmic [26] produced analytic expressions for high SNR 

1 12 (wT s ) 4 

a/> — — 

1 T s 2 M [\-\2 (wTS) 2 ] 


(37) 


and low SNR 

1 4Vtc(h'T 5 ) 3 w 2 
T s 2 M S 2 



(38) 


B. Cramer-Rao Bound for Mean Frequency Estimation 

Many different estimators for mean-frequency have been suggested and analyzed. In the limit as 
M—>°» the performance of the maximum likelihood (ML) estimate approaches the CRB [10,20,31,32] and 
the estimates have a Gaussian distribution. However, the ML estimate is a numerically intensive algo- 
rithm. The performance of an approximate ML estimator will be compared to the exact and approximate 
CRB’s. To reduce the computational effort, the spectral width w is assumed known and the signal power 
S is estimated by 

S = \-^ 7 '£[P(m)-N]\ (39) 

M m=0 


and | . | denotes absolute value. An approximate maximum likelihood estimate for the mean frequency 
that assumes knowldege of the spectral width w (this is true when the velocity fluctuations over the sens- 
ing volume of the pulse are negligible and the spectral width is determined by the transmitter pulse) is the 
value / [10] that maximizes the log likelihood function L assuming the periodogram elements are mutu- 
ally uncorrelated, i.e., 


u y l PQn) 
hP(m t Sf) ‘ 


(40) 


Note that the average periodogram P(m) is used instead of the continuous spectrum [10,17]. This 
improves performance for narrow spectral width and high SNR. A version of Levin’s ML estimator that 
uses accumulated periodograms has been investigated^ by Rye and Hardesty [17,18]. To produce a 
numerically efficient estimator, L is calculated for / = k f Af ,kf =1,2,. M . The performance of this 
approximate ML estimator was determined from 10000 simulations of random data [3] with M =128 
using the covariance model Eq. (29) with a mean frequency f=Q.5/T s , i.e., / is located in the center of 
the periodogram bin k=MI2. The Probability Density Function (PDF) of the mean frequency index kf of 
these approximate ML estimates is shown in Fig. 1 for SAT? =0.1. The PDF’s are characterized by a uni- 
form distribution of "bad" estimates and a clump of "good” estimates about the true mean frequency 
[3,4,1 1,12,17]. Because the signal power for each estimate is random, there will be realizations where the 
signal will be hidden in the additive noise and a random noise spike will be chosen as the mean fre- 
quency. The "bad" estimates produce a bias which creates a problem for comparison with the CRB calcu- 
lated for unbiased estimates. One solution is to calculate die CRB for a biased estimator [20,31], but the 
bias depends on the estimator. This comparison is presented in [3]. We will compare the CRB to the 
standard deviation of the good estimates around the true mean frequency. A theoretical PDF is calculated 
assuming two components: a uniform distribution for die "bad” estimates and a Gaussian distribution for 
the "good” estimates. The standard deviation of the good estimates and the fraction of bad estimates are 
determined by minimizing the mean square error between the theoretical PDF and the simulated PDF. 
For the parameters of Fig. 1. and a normalized spectral width wTgMkOl, 25.1% of the estimates are bad 
estimates; for wTjsO.l, 28.8% of the estimates are bad estimates. 

A comparison of exact and approximate CRB’s for the estimation of the mean frequency is shown 
in Fig. 2 for A/=128 data points, as a function of the normalized spectral width wT s for various SNR. 
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When the maximum correlation p miX is small (p mJ o t <0.01), the periodogram CRB is a good approxima- 
tion to the exact CRB. When the spectral width w is less than the periodogram resolution A / , the 
periodogram CRB is larger than the exact CRB and p ra „ approaches unity. This is the explanation for the 
failure of the periodogram CRB and all it’s analytic approximations (Whittle’s, Levin’s, Zmic’s) to agree 
with the exact CRB in the limit of spectral width approaching zero, (a sine wave with random amplitude 
and phase) (see Ref [26] Eq. (A. 18)). Levin’s approximate CRB is slightly lower than the exact CRB for 
narrow spectral width and much lower at large spectral width. 

The performance of the approximate ML estimator is complicated at low SNR by the bias intro- 
duced by the bad estimates (see Fig. 1). For SNR=0.1, the standard deviation of the good estimates is 
below the exact CRB for large spectral width. However, the standard deviation using all the estimates is 
above the exact CRB. For higher SNR, there is little bias because most of the estimates are distributed 
around the true mean frequency and the standard deviation of the approximate ML estimates approach the 
exact CRB for moderate to large spectral width (see also Ref. [1]). For very narrow spectral width all the 
signal power is located in one periodogram element if e=0. If the true mean frequency is not centered on a 
periodogram element (e*0), then the signal power is distributed over more than one periodogram element 
[2, Chap. 5] and the approximate ML estimate should include the dependence on e. This regime is rarely 
encountered in Doppler radar but is important for Doppler Lidar measurements when the wind fields over 
constant the sensing volume of the pulse. Then the exact maximum likelihood estimator is preferred [3]. 

A comparison of the exact and approximate CRB’s for the estimation of the mean frequency / are 
shown in Fig. 3 versus SNR for M =128 data points and various spectral widths. Again, the approximate 
periodogram CRB agrees with the exact CRB when the maximum correlation p max <0.01. Levin’s approx- 
imate CRB is a good approximation to the exact CRB for this parameter regime. Zmic’s low-SNR ana- 
lytic expressions for CRB (a) are accurate for only very low SNR. The high-SNR analytic approximation 
of Zmic (b) [Eq. (37)] does not agree with the exact CRB due to the poor representation of the summation 
by integration (see Fig. 4). The standard deviation of the approximate ML estimates is in excellent agree- 
ment with the periodogram CRB. This result appears to be general [3], i.e., the standard deviation of 
mean frequency estimates based on the periodogram or based on estimates of the covariance are bounded 
by the periodogram CRB if there is no bias from the bad estimates (higher SNR). This class of mean fre- 
quency estimators includes the periodogram estimators, auto-regressive estimator (AR), auto-regressive 
moving average (ARMA), and Minimum Variance or Capon spectral estimator [7]. 

The source of the error in the high SNR analytic expressions for CRB and Levin’s approximation 
are revealed in Fig. 4 for a very high SNR of 1000. The terms of the Fisher information matrix element 
J'n with 6=0 [see Eqs. (16), (26), (34), (35)] are compared with Levin’s approximation (see Eq. (A. 36) of 
Ref. [26]) and the high-SNR approximation of Zmic [26]. The high SNR approximation for the integrand 
is dominated by the largest values of the periodogram index m , i.e., the periodogram coefficients farthest 
removed from the spectral peak of P(m). This region is the most difficult to represent with analytic 
approximations to the discrete summation, especially for small spectral width. This error is not reduced as 
SNR approaches infinity. Levin’s approximation to Eq. (26) is also higher and this produces a CRB that 
is lower than the periodogram CRB but this error results in a good approximation to the exact CRB for 
much of the parameter space of interest. 


C. Cramer-Rao Bound for Signal Power and Spectral Width 

Estimators for the signal power and spectral width have been discussed previously ([2], [9], [10], 
[13], [26], [27]). The exact CRB for estimation of the signal power when spectral width is also estimated 
is given by 


o$> 


J 33 


J 2hf 33 ~J%3 


(41) 


The CRB’s for estimation of S is shown in Fig. 5 versus normalized spectral width for various SNR. The 
periodogram CRB is approximately equal to the exact CRB when the maximum periodogram correlation 
Pm«x > s small. Levin’s approximate CRB is poor for small and large spectral width. For low SNR, the 
CRB is close to the value unity which is expected for one independent sample of the signal per estimate. 
For high SNR, the signal power can be more accurately estimated when the spectral width is large. This 
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reflects the larger number of independent samples available per estimate. 

The CRB’s arc shown in Fig. 6 versus SNR for small and moderate spectral width. The exact CRB 
decreases with increasing SNR. The periodogram CRB is a poor approximation for the exact CRB for 
narrow spectral width and high SNR. This is reflected in the high correlation p m „. The approximate ana- 
lytic calculations agree at low SNR (a) but not at high SNR (b). 

The exact CRB for estimates of the spectral width when signal power is also estimated is given by 




J 11 

J 11-^33 ?3 


(42) 


The CRB’s for estimates of spectral width w are shown in Fig. 7 versus normalized spectral width wT s 
for various SNR. The ability to estimate spectral width using the periodogram degrades when the spectral 
width is smaller than the frequency resolution tif (wMT s <l). The periodogram CRB agrees with the 
exact CRB when the maximum correlation p m4X is small. Levin’s approximate CRB is less than the exact 
CRB for small spectral width. 

The CRB’s for estimates of w are shown in Fig. 8 versus SNR for small and moderate spectral 
width. The exact CRB decreases with increasing SNR. Again, Zmic’s low SNR limit (a) agrees with the 
exact CRB at low SNR but Zmic’s high SNR limit (b) is incorrect The source of the error at high SNR is 
the failure of the integral approximation for the discrete summation (see Fig. 4). The approximate CRB’s 
for signal power and spectral width are a good approximation to the exact CRB’s when pmax<0.01 for all 
cases considered. 


IV. SUMMARY 

Previous results for the exact CRB of the parameters of a complex Gaussian random process are 
applied to Doppler radar/lidar parameter estimation. A new expression Eq. (15) for the exact CRB for 
estimation of mean frequency was derived for a complex Gaussian signal in uncorrelated noise. For this 
model, the CRB’s for all the signal parameters are independent of the mean frequency. Approximate 
CRB’s were derived by taking advantage of the fact that the DFT of the data has the same CRB and there- 
fore contains the same information as the original time series. When the DFT coefficients are mutually 
uncorrelated then so are the periodogram coefficients, and the periodogram CRB is equal to the exact 
CRB. For this case, the periodogram (or equivalently, the estimates of the covariance) contain the same 
information as the original time series. The periodogram CRB and Levin’s approximate CRB were com- 
pared to the exact CRB for estimation of mean frequency, signal power, and spectral width for a signal 
with a Gaussian covariance function. Levin’s approximate CRB is lower than the exact CRB for narrow 
spectral width. For moderate spectral width and moderate SNR, all the CRB’s show good agreement 
The performance of an approximate maximum likehood estimate for mean frequency / Eq. (40) using the 
periodogram and assuming knowledge of the spectral width approaches the exact CRB for moderate SNR 
and moderate spectral width. The role of the CRB for low SNR is complicated because all estimates for 
mean frequency are biased due to the uniformly distributed "bad" estimates (See Fig. 1 and Ref. [3]). The 
periodogram CRB appears to be a lower bound for the unbiased spectral domain estimates of the mean 
frequency [3] (periodogram, maximum entropy, auto-regressive, auto-regressive moving average, 
m inimum variance spectral estimator (Capon) [7]). For narrow spectral width and high SNR, more than a 
factor of two improvement in mean-frequency accuracy is indicated by the CRB [see Fig. 2]. This 
impro vement can be realized with the Maximum Likelihood estimator [3]. The calculation of the 
periodogram CRB and the maximum correlation coefficient p™** can be efficiently performed using the 
FFT algorithm. Previous analytic high SNR limits for exact CRB’s were shown to be inaccurate because 
the discrete summations cannot be accurately approximated by integration (see Fig. 4). 
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FIGURE CAPTIONS 


Fig. 1 Probability Density Functions (PDF) (•) from 10000 approximate Maximum Likelihood estimates 
of mean frequency / Eq. (40) for a SNR =0. 1 and two different spectral widths. The best fit 
theoretical model is indicated by ( ). 

Fig. 2 CRB for mean frequency / vs normalized spectral width for SNR of 0.1,1, and 10 and Af =128 
data points. Exact CRB [Eqs. (16), (1 1), and (30)] ( ), periodogram CRB (....), Levin’s CRB (- 

- -), and maximum periodogram correlation p m „ (._._). The performance of the approximate 
maximum likelihood estimator Eq. (40) is indicated by O. 

Fig. 3 CRB for mean frequency / vs SNR for normalized spectral widths of 0.01 and 0.05 and M=128 
data points. Exact CRB [Eqs. (16), (1 1), and (30)] ( _J, periodogram CRB (....), Levin’s CRB (- 

- -), Zmic’s low SNR limit Eq. (38) a, Zmic’s high SNR limit Eq. (37) b, and maximum periodo- 
gram correlation p max (._._). The performance of the approximate maximum likelihood estima- 
tor Eq. (40) is indicated by O. 

Fig. 4 Integrand for approximate CRB’s for mean frequency / vs periodogram index m for an 
SNR=1000, M =128 data points, and normalized spectral widths of 0.1 and 0.05. The elements of 
the discrete summation for the Fisher information matrix element J' n Eq. (26) with Eqs. (34), 

(35), and m F =0 are indicated by ( ), Zmic’s analytic approximation for Levin’s CRB [Ref. 

[26], Eq. (A.36) is (....), Zmic’s high SNR approximation [26] is ( ), and the average periodo- 

gram Eq. (34) is marked as P (m ). 

Fig. 5 The exact CRB for signal power [Eqs. (41), (11), and (30)] ( ), the periodogram CRB (....), 

Levin’s CRB ( ), and the maximum periodogram correlation p max (._._) vs normalized spec- 

tral width wTs . 

Fig. 6 The exact CRB for signal power [Eqs. (41), (11), and (30)] ( ), the periodogram CRB (....), 

and Levin’s CRB ( — ), vs SNR for normalized spectral widths of 0.01 and 0.05 and M=128 data 
points. Zmic’s [ Ref. [26], Eq. (A.40)] low SNR limit a, and Zmic’s [ Ref. [26], Eq. (A.38)] high 
SNR limit b. The maximum periodogram correlation p m „ is (. _ . _). 

Fig. 7 The exact CRB for spectral width [Eqs. (42), (11), and (30)] ( ), the periodogram CRB (....), 

Levin’s CRB ( — ), and the maximum periodogram correlation p mM (._._) vs normalized spec- 
tral width wT s for various SNR and M =128 data points. 

Fig. 8 The exact CRB for spectral width [Eqs. (42), (11), and (30)] ( ), he periodogram CRB (....), 

Levin’s CRB ( — ), and the maximum periodogram correlation p m „ (._._) vs SNR for normal- 
ized spectral widths of 0.01 and 0.05 and M=128 data points. Zrnic’s [Ref. [26], Eq. (A.41)] low 
SNR limit a, and Zmic’s [Ref. [26], Eq. (A.39)] high SNR limit b. 
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Abstract 

The performance of a Coherent Laser Radar is determined by the statistics of the 
coherent Doppler signal. The heterodyne efficiency is a excellent indication of per- 
formance because it is an absolute measure of beam alignment and is independent 
of the transmitter power, the target backseat ter coefficient, the atmospheric attenu- 
ation, and the detector quantum efficiency and gain. The theoretical calculation of 
heterodyne efficiency for an optimal monostatic lidar with a circular aperture and 
Gaussian transmit laser is presented including beam misalignment in the far-field and 
near-field regimes. The statistical behavior of estimates of the heterodyne efficiency 
using a calibration hard target are considered. For space based applications, a biased 
estimate of heterodyne efficiency is proposed that removes the variability due to the 
random surface return but retains the sensitivity to misalignment. Physical insight 
is provided by simulation of the fields on the detector surface. The required detector 
calibration is also discussed. 

Key words: Lidar, Laser Radar, Remote Sensing. 


1. Introduction 

An important measure of the performance of coherent laser radar (CLR) is the Signal-to- 
Noise Ratio (SNR) defined by [ 1] 

SNR(t) = (1) 

(*%) 

where () denotes ensemble average, (i|(t))[i4 2 ] is the average coherent power at time t[s], 
is the average noise power. For an infinite diffuse target at range /?[m] with 
uniform backscatter coefficient p, detector noise dominated by Local Oscillator (LO) shot 
noise, a detector with uniform quantum efficiency tjq that collects all the energy of the LO 
beam, the CLR SNR is given by [ 1, 2] 

SNR(R) = * Q fy~ C{R) (2) 

where Pl\W] is the transmitted laser power, h[Js] is Plancks constant, is the frequency 

of the laser, 5[s -1 ] is the signal bandwidth of the detector, K is the one-way atmospheric 
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attenuation, and C(R ) is the coherent responsivity. For an aerosol target, similar expres- 
sions are obtained [Ref. [ 1], Eqs. (89), (92)]. When the effects of atmospheric refractive 
turbulence are negligible 

C(R ) = ArTtT}h(R)/R 2 (3) 

where v4/*[m 2 ] is the area of the telescope aperture, TV is the fraction of the laser power 
transmitted through the telescope aperture, and tjh{R) is the heterodyne efficiency. Het- 
erodyne efficiency is unity when the backscattered field is proportional to the LO field on 
the detector surface. Heterodyne efficiency is a useful measure of performance for a CLR 
because it is independent of the transmitter power, the backscatter coefficient, the detector 
quantum efficiency, and the atmospheric attenuation. The heterodyne efficiency can also be 
estimated from the detector signal [ 1], i.e., 

jjlOO) (4) 

2(lDc)(Is(t)) 

where Idc[A ] is the direct current (DC) due to the LO, i s intermediate frequency 

(IF) signal current or coherent signal, and 7s(<)[A] is the direct detection signal current due 
to the backscattered field [Note the typo in Ref. [ 1], Eq. (15)]. For diffuse and aerosol 
targets, the backscattered field at the receiver is spatially random and described by "speckle” 
statistics, i.e., the complex field has a zero- mean joint Gaussian probability density function 
(PDF). The statistical behavior of estimates for the heterodyne efficiency will be investigated 
with theoretical calculations and numerical simulations. 



2. Theory 

The analysis of heterodyne efficiency is presented using the normalized direct detection 
power 

d = I s I (Is) „ (5) 

and the normalized coherent detection power 

C = 2(/Jm 7s) (6) 

Then (d) = 1 and (c) = rjff. 

The heterodyne efficiency depends on the CLR design, which can be represented by [ 
1] the transmitter aperture function WV(u), the normalized transmitter laser field incident 
on the transmitter aperture e£,(u, 0)[m _1 ], the receiver aperture function HVj(v), and the 
normalized Local Oscillator (LO) field at the receiver aperture e£,o(v, 0)[m -1 ]. Here, u[m] 
and v[m] denote the two-dimensional transverse vector for the transmitter aperture and 
receiver aperture, respectively. 

We will consider a monostatic CLR [WV(u) = WVi(u) = W(u)] with a circular telescope 
aperture of radius a[m], i.e., 


2 



W{ u) = 1 |u| < a 

= 0 |u| > a 

and 

An = f |iy(v)| 2 dv = 7r a 2 = 7rD 2 /4 

where Z)[m] is the diameter of the aperture. 

The heterodyne efficiency can be written as [ 1] 

/?2 \2 >oo 

Vh{R) = R)jBPLo(P, R)dp 

where 

ir(p,ii) = |er(p, f?)| 2 


(7) 


( 8 ) 


( 9 ) 


( 10 ) 


jBPLC)(P,R ) = |CBPL0(P, 7Z)| 2 


( 11 ) 


are the irradiance [m 2 ] of the normalized transmitter and Back Propagated Local Oscillator 
(BPLO) field, respectively. Here, 


e r (p ,R) = f ex(u, 0)G(p; u, R)du 

is the normalized transmitter field at the target coordinate (p, R) [m], 
e T (u,0) = e L (u,0)fT(u) 

is the normalized transmitter field at the telescope aperture, 
e BPLo(p, R) = f cbplo(v, 0)G(p; v, R)dv 
is the normalized BPLO field at the target, 
c-BPLoiy , 0) = e2(v,0)fT(v) 

is the normalized BPLO field at the receiver aperture, 


G(p;u,R) = 


2niR 


exp 


z k 


( 12 ) 


(13) 


(14) 


(15) 


(16) 


is the free space Green’s function [ m ~ 2 ] for propagating the scalar fields from the transmitter- 
receiver plane to the target plane under the Fresnel approximation where fc[m -1 ] = 2x/A is 
the wavenumber of the transmitter field, and A[m] is the wavelength of the field. 

The integrated normalized irradiance 


/ °° 

Jr(p, R)dp = Tj 

•OO 


(17) 
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is equal to the power transmittance T t of the laser through the transmitter-receiver aperture. 
Similarly, 

f°° j B PLo(p,R)dp = T r (18) 

J—oo 

is equal to the power transmittance Tr of the BPLO field through the transmitter-receiver 
aperture. 

A Gaussian spatial distribution is the most common laser transmitter field, i.e., 

e t (u,0) = (2™i)- l/, exp (-^j - ^0 (19) 

where a L [m] is the 1/e intensity radius of the beam at the transmitter aperture, and Fi,[m] is 
the radius of curvature of the phase front. The irradiance at the target is given by [substitute 
Eqs. (13) and (19) into (12), convert to polar coordinates and integrate over angle] 

= < 2 °) 

■KO h 

where 

ft* = 6 2 - iS 2 RF ( 1 - R/F l ) (21) 

St = a/vL > &RF = o/Rf , Rf = (R/k) 1 ^ 2 (22) 


Y(Cl,q) = / exp(— x 2 (l 2 /2)Jo(qx)xdx 
Jo 

An efficient series expansion for Y(Cl,q) is [ 3] 


r(ft,g) = n -2 

[ exp ( 2fl ! ) exp ( T^f-G?) Jn(,) ] 

= ft -2 1 

00 /ft 2 \ n 

: xp(-n J /2) 2 (— j jm 


q < |0 2 | 
q > |ft 2 | 


(23) 


(24) 


The far- field result is produced when ftj = S F , which occurs when R = Fl (focused condi- 
tion). For a collimated beam (Fl = oo), the far-field condition is Sj » Srf or Rf » <?l- 
A Gaussian field focused on the detector is a common choice for the LO field. The 
receiver- plane representation of the LO field is 

eLo(v, 0) = (2 tt<7 lo)~ 1/2 exp (“7^ + ( 25 ) 


and 

]bpw(v,R) = W *)! 2 


(26) 
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where 


= S R — iS RF (l — R/Flo) , Sr = a/c^LO (27) 

In practice, it is difficult to perfectly align the LO field with the backscattered field. The 
misalignment angle A 0{rad\ produces an offset A p — RAd[m\ between the two irradiance 
distributions in the target plane. (The effects of refractive turbulence on beam misalignment 
for Gaussian beams and a Gaussian aperture were investigated by Frehlich [4]). For this 
CLR, the heterodyne efficiency including misalignment is given by 

m (R) = i|S<3(n T ,n s ,A,) (m) 

where 

/ OO 

|K(nr,q)| J |K(n«,q-Aq)|Vq (29) 

-CO 

is the overlap integral of the normalized transmitter irradiance |y(flr,q)| 2 and BPLO irra- 
diance |K(fIn,q — Aq)| 2 and 

A q = ttDAO/X (30) 


is the normalized angle-misalignment. 

Wang [ 5] calculated the optimal parameters for a Gaussian transmit beam and a Gaus- 
sian LO beam focused on the detector for the far-field condition [Sj — Sr — 1.7633, 
t)h = 0.41992, Tt = Tr = 1 — exp(— 6\) = 0.95536). Rye [ 6] proposed an optimal de- 
sign where the BPLO beam does not have any truncation by the telescope aperture, i.e., 
Tr = 1. For this design 

w(fl)= -I, 4 ***” , 2 ,. Q((tT,fl R ,Aq) (31) 

7rr r [l - exp (—<%)] 

and for a Gaussian LO beam, the optimal parameters axe St = 1.7367, Tj = 1 —exp (—Sj) = 
0.95101, Sr = 1.1931, t)h — 0.460958. The difference between these two designs is discussed 
by Rye and Frehlich [ 7). 

For a collimated lidar, the near field behavior of heterodyne efficiency can be approxi- 
mated by 

\2 yoo 

Vhnf(R) = ^ R j> T J_ jr(p,0)j B PLo(p,0)dp (32) 

For the collimated Gaussian lidar with the Rye design [see Eqs. (11), (13), (15), (16), (19), 
(25), (32) and normalize the BPLO field to produce Tr = 1] 


(R , 4ffiw[l - exp(— - SQ j 

Vhnf{ 8 4 RF (Sl + <%)[! - exp(-^)][l - exp(-<%)] 


(33) 


The optimal LO field for a Gaussian transmitter field and a circular aperture was deter- 
mined by Frehlich [ 8] and produces essentially the same performance as the design proposed 
by Rye. The behavior of heterodyne efficiency for a CLR with the Rye design will be con- 
sidered here. 
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3. Simulation of Performance 


One important aspect of CLR performance for diffuse and aerosol targets is the random 
nature of the backscattered field. The superposition of many random backscattered fields 
with random amplitude and phase produces a complex Gaussian distribution or "speckle 
field” for the total backscattered field [ 9]. The detection of such random signals in random 
noise has been considered by Shapiro et al. [ 10]. We will investigate the estimation of 
the heterodyne efficiency by generating computer simulations of the normalized heterodyne 
power c and the normalized direct detection power d. Realization of the fields on the detector 
will also be calculated to provide insight. 

The backscattered field at the target can be represented as 


e s (p,R) = U(p)e T (p,R) (34) 

where U( p) is a random complex amplitude with the property 

(C/(p 1 )f/*(p 2 )) = A 2 ^(p 1 -p 2 ) (35) 

where £(p)[m -2 ] is the two-dimensional delta function. The backscattered field at the re- 
ceiver plane es(v, 0) is given by 

es(v,0) = / e s (p,R)G{p;v,R)dp (36) 

The normalized direct detection power is then 

d = ^ £V(v)| J |e s (v,0)| 1 <lv (37) 

and the normalized coherent detection power is 

y*00 

C = ~A^Ft l L es{ *' ( 38 ) 


The effects of angle misalignment A 6 are calculated by multiplying the BPLO field with 
no misalignment by the phase term [ 4] exp(ikA6v x ) where t> r [m] is the x-coordinate of 
the vector v. The random fields are calculated using the Fast Fourier Transform (FFT) 
algorithm to propagate the fields. The random amplitudes U(p) are generated using the 
Box- Muller algorithm. 

A typical receiver employs optical elements that pass all the backscattered field. Then, 
for ideal optical elements, the field on the detector surface is 

e s (w, F) = ^^(vjesiv, 0) exp(- jv • w )dv (39) 

where F[m\ is the effective focal distance of the optical system. The fields on the detector sur- 
face are particulary intuitive since the heterodyne efficiency is unity when the backscattered 
field is proportional to the LO field. This is also true for the receiver-plane representation. 
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4. Results 


The simulation results for the far-held optimal design of Rye and the corresponding col- 
limated near-held case with 6rf = 4 will be presented. The statistical variability of the 
coherent detection power and the direct detection power from a simulation of 10,000 shots 
for the perfectly aligned CLR for the far-held case is shown in Fig. 1 and for the near-held 
case in Fig. 2. If the heterodyne efficiency had been unity, all the points would lie on a 
straight line with a slope of unity. The far-held case has more scatter than the near held 
case. This is due to the increased spatial averaging of the backscattered irradiance by the 
telescope aperture in the near-held case. 

Typical realizations of the backscattered held on the detector surface for a perfectly 
aligned CLR is shown in Fig. 3 for the far-held and near-held cases. The dimensions of 
the central lobe is approximately the diffraction limit A F/a of the circular aperture. The 
poor match between the backscattered held and LO held for the typical far-held case is 
the reason the heterodyne efficiency is 0.46096 for the perfectly aligned far-held case. More 
than half the coherent power is lost due to the random mismatch between the backscattered 
and LO helds. The near-held case is characterized by a larger distribution of the held 
over the detector surface with huctuations on the order of the diffraction limited scale but 
with reduced amplitude which results in a heterodyne efficiency of 0.021832. The ratio of the 
heterodyne efficiency to the approximate near-held heterodyne efficiency Eq. (33) is 1.05554. 
The large number of "speckles” or regions of coherent held produces smaller variability of 
the direct detection signal (see Fig. 2). Also, for typical CLR designs, the detector would be 
smaller than the spatial extent of the held and the direct detection power would be reduced. 
Measurements of heterodyne efficiency in the near-held regime may contain a bias due to 
the finite extent of the detector. A rare case of high heterodyne efficiency and very low 
heterodyne efficiency is shown in Fig. 4 for the far-held case. 

The sensitivity of heterodyne efficiency to the misalignment angle Ad of the transmitter 
and BPLO axis [ 1, 6] is shown in Fig. 5 for the far-held and near-held cases. The theoretical 
heterodyne efficiency depends on the single parameter A q = kaAd = nDAO/ A. The simu- 
lation results are in excellent agreement with the theory for both the near-held and far-held 
case. (The estimation error for each simulation data point is 1%, consistent with exponential 
statistics for the coherent power.) The heterodyne efficiency was estimated by the ratio of 
the average of c to the average of d using the 10,000 shots of simulated data in Figs. 1 
and 2 and various angle misalignment. In the laboratory one would measure the heterodyne 
efficiency in this way using a calibration target with a uniform backscatter coefficient p. The 
statistical accuracy of the estimate for heterodyne efficiency is determined by the number 
of shots averaged and the statistical distribution of the two estimates for < c >= tj h and 
< d >= 1. The Probability Density Function (PDF) of the direct detection power d and 
the Gamma Distribution function g(d] a , /?) proposed by Goodman [ 9] are shown in Fig. 6. 
The Gamma distribution function is given by 

9{d-,a,/3 ) = ^-pjd 0 - 1 exp{-f3d) (40) 

and the parameters of the Gamma Distribution are a = /? = where <r$ is the variance of 
d. For the far-held case a* — 0.531. The same comparison for the near- held case is shown in 
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Fig. 7. Here, 04 = 0.1526. In both cases the simulation results are in reasonable agreement 
with the Gamma distribution. 

The PDF of the normalized coherent detection power for the far-field case is shown in 
Fig. 8 for perfect alignment along with the theoretical exponential distribution. Because the 
backscattered field is the superposition of many, randomly phased components, the coherent 
power should have an Exponential PDF. This was verified by the simulations for all cases. 

The estimate of heterodyne efficiency r\u with a fixed calibration target using the average 
of J shots is given by 


eL, 


EL, d„ 


N_ 

M 


(41) 


where c* and dk are the coherent power and direct detection power for shot k. The numerator 
N has a Gamma Distribution g(N;l/(c),J) and the denominator M also has a Gamma 
Distribution g(M;l/a%, J/crl). For large J, the statistical accuracy of the estimate rju is 
approximately 


A f)n _ AN AM _ 1 + <Td 
~^T ~ ~N~ + ~M~ ~ v7 


(42) 


The estimation of heterodyne efficiency for satellite based lidars is a promising technique 
to verify and maintain beam alignment, but is more difficult since there is no calibration 
target. The surface return is the most promising signal since it may be large enough to 
have a detectable direct detection signal but there is added variability from shot to shot 
due to the random backscatter coefficient of the earths surface and the random atmospheric 
attenuation. This added variability produces larger scatter in the estimate of the coherent 
and direct-detection signed power. However, for estimation of heterodyne efficiency, this 
variability can be removed by defining a biased estimate of heterodyne efficiency for each 
shot as 


b = c/d (43) 

provided the direct detection signal can be estimated accurately for every shot. The average 
of this biased estimate and it’s standard deviation as a function of misalignment angle A 9 
is shown in Fig. 9 along with the true heterodyne efficiency (see Fig. 5). The bias of the 
average of b is small and therefore provides a useful absolute measure of beam alignment. 
The statistical reliability of this biased estimate depends on the PDF of 6, which is shown 
in Fig. 10 for perfect alignment. Also shown is the PDF of the average of 10 estimates 
of b and a Gaussian PDF with the correct mean value and variance given by the variance 
of 6 divided by the number of shots averaged. The same results for a misalignment of 
A q = kaAO = ttDAO/X = 2 which corresponds to 2.08 dB loss in SNR are shown in Fig. 
11. The Gaussian approximation for the PDF of the average of 10 estimates of b is a good 
approximation. However, am accurate estimate of the biased heterodyne efficiency requires 
an average of many shots. The PDF of a 40 shot average of b is shown in Fig. 12 for perfect 
alignment and misalignment of irDA9/\ = 2 (2.08 dB loss in SNR) and irDAO/X = 3 (4.76 
dB loss in SNR), using the Gaussian approximation for the PDF’s. The CLR misalignment 
can be reliably determined with 40 surface returns. 


8 



5. Detector Calibration 


A real doppler lidar system does not measure the signal currents directly. Therefore, detector 
calibration is required to extract the required signal currents from the detector signal. The 
heterodyne efficiency can be estimated using both Continuous Wave (CW) or pulsed mode 
of operation. For CW operation, two signal currents Is(t) and Ipc are DC signals which 
are given in terms of the detector voltage V by 

I(V) = G D h(V) (44) 

where /(V) is the photodetector response function G^A/V] is the constant for a linear 
response, and A( V r )[V r ] is the nonlinearity function or deviation from linearity. We assume 
that the photodetector is ideal in the limit of small voltage, i.e., h(V) = V. For pulsed 
operation, the frequency response of the detector may be required to estimate the direct 
detection signal /s(f). 

If the magnitude of the IF signal current is(t) is small enough and if the magnitude of 
the IF voltage signal t>s(f)[K] does not change with the IF frequency v[Hz] then 


*s(0 = 

G D f(V)v s (t) 


(45) 

where 




m = 

1 dI(V) 
G d dV 

dh(V) 

dV 

(46) 


is proportional to the derivative of the photodetector response function and /( 0) = 1. A 
measurement of heterodyne efficiency requires a measurement of h(V) and f(V), which can 
be performed by the method proposed by Frehlich [ 11]. (Note that Ref. [ 11] considered 
the calibration of integrated irradiance over the detector surface which is proportional to the 
detector current). The operating point (7 0 , Vo) determines the relationship between detector 
voltages and the ideal signal currents. The direct detection currents are given by 

Idc = Goh(\ o) (47) 

and 

Is(t) = G D f(V 0 )Vs(t) (48) 

where Vs(f)[V] is the increase in the DC detector voltage due to the backscattered field 
collected by the detector. 

The effects of the detector electronics and detector quantum efficiency as a function of 
IF frequency are described by the normalized frequency response function 


H(-) = 

ws(» / ) 

(49) 

vs(v — *• 0) 

Then 

*5(0 = 

G D f(V Q )v s (t)/H(i v) 

(50) 
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and the heterodyne efficiency is given by [see Eq. (4)] 


/WfrKO) 

2h(V 0 )(V s )H*(v) 

which is expressed in terms of the calibration functions /(Vo), h(Vo), and the detector 

voltages i>s(f) and Vs(t). Note that the ideal detector gain Go is not required. The sensitivity 
of the calibration functions to the detector operating point (7o, Vo) should be investigated. 




0. Summary 

The estimation of heterodyne efficiency is an excellent measure of CLR performance since 
it is proportional to SNR but independent of the transmitter power, the target backscatter 
coefficient, the atmospheric attenuation, and the detector quantum efficiency and gain. The 
heterodyne efficiency is an absolute measure of the alignment of the backscattered and LO 
field on the detector surface. The estimation of heterodyne efficiency is a statistical problem 
because the detector signals are random. An estimation of the heterodyne efficiency requires 
an accurate average of the coherent detection power, the direct detection power, and the LO 
current. When using a calibration target, an unbiased estimate of the heterodyne efficiency 
can be obtained by an ensemble average of the coherent and direct detection power. For a 
lidar in space, the ratio of the coherent detection power to the direct detection power for each 
ground return is a convenient biased estimate of heterodyne efficiency. The average of this 
biased estimate produces a sensitive measure of CLR alignment and removes the statistical 
variability of the random fluctuations of the signals due to the random backscatter coefficient 
of the earth surface and the random atmospheric attenuation. Using 40 shots produces 
an effective measure of system performance and beam alignment [see Fig. 12]. All the 
results presented here assume a large detector that collects all the LO power and all the 
backscattered power. 

The calibration of real detectors requires two measurements: the low frequency response 
function h(V) (which is obtained [ 11] from a measurement of /(U), the derivative of h(V)) 
and the frequency response H{v). 
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Coherent Detection Power c 


Figure 1 

Fig. 1. Normalized coherent detection power c versus the normalized direct detection power d 
for a perfectly aligned CLR and 10,000 shots with the target in the far field. 
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Figure 2 

Fig. 2. Normalized coherent detection power c versus the normalized direct detection power d 
for a perfectly aligned CLR and 10,000 shots with the target in the near field (6 RF = a/Rr = 4). 
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Fifi- Simulation of the backseat tered field on the detector surface in the x direction F) 

and the y direction e s (w y , F) for typical far field ( t) H = 0.460, d = 1.015, c = 0.467) and near field 
with 6 ftp = a/Rp = 4 ( 77 // = 0.0785, d = 0.969, c = 0.0760) cases. The real part of the complex 
field is ( — ) and the imaginary part is (• • •). The optimal LO field e LO {w , F) is also shown. 
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Fig. 4. Simulation of the backscattered field on the detector surface in the x direction es(u? x , F) 
and the y direction es(w y , F) for far field and rare cases of excellent heterodyne efficiency 
= 0.9465, d = 3.325, c = 3.147) and poor heterodyne efficiency (t]h = 0.000344, d = 0.830, 
c = 0.000286). The real part of the complex field is ( — ) and the imaginary part is (• • •). The 
optimal LO-field eio{w,F) is also shown. 





Heterodyne Efficiency r\ 


Figure 5 


Fie 5 Heterodyne efficiency W as a function of misalignment *D A0/A with the results 
of ihe 10,000 shot simulation of Fig. 1 for the far-field (•) and near-held (o) cases with 
6 rf = a/R F = 4. The theoretical calculations Eq. (31) are ( — ). The best-fit Gaussiaui functm^ 
0.46096 exp[-(fl’.DA0/2.8363A) 2 ] is (• • •) and the best fit power-law model O.9157(7r£>A0/A) 

is (---)• 
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Figure 6 


Fig. 6. PDF of direct detection power d from simulation of Fig. 1 and the theoretical Gamma 
distribution (— ) for a perfectly aligned CLR in the far-field. 
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Figure 7 

Fig. 7. PDF of direct detection power < 
1 and the theoretical Gamma distributio 
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Figure 8 

Fig. 8. PDF of coherent detection power c from simulation of Fig. 1 and the theoretical Expo- 
nential distribution ( — ) for a perfectly aligned CLR in the far-field. 
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Fig. 11. PDF of the biased estimate of heterodyne efficiency b = c/d for a misalignment of 
ttDAO/X = 2 in the far-field for 1 shot and the average of 10 shots. The Gaussian PDF ( — ) with 
the estimated mean < b > and standard deviation <Jt>. 
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Figure 12 


Fig. 12. PDF of the average of 40 estimates of b for a misalignment of irDAd/X — 0 ( ), 2 

(• • •), and 3 ( ) assuming a Gaussian PDF. 
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Abstract 

The performance of coherent Doppler lidar is determined by the statistics of the 
coherent Doppler signal. The derivation and calculation of the covariance of the 
Doppler lidar signal is presented for random atmospheric wind fields and wind shear. 
The signal parameters are defined for a general coherent Doppler lidar system in 
terms of the atmospheric parameters. There are two distinct physical regimes: one 
where the transmitted pulse determines the signal statistics and the other where 
the wind field and atmospheric parameters dominate the signal statistics. When the 
wind fields dominate the signal statistics, Doppler lidar data is non-stationary and the 
signal correlation time is proportional to the operating wavelength of the lidar. The 
signal covariance is derived for single shot and multiple shot conditions. For a single 
shot, the parameters of the signal covariance depend on the random, instantaneous, 
atmospheric parameters. For multiple shots, various levels of ensemble averaging over 
the temporal scales of the atmospheric processes are required. The wind turbulence 
is described by a Kolmogorov spectrum with an outer-scale of turbulence. The effects 
of the wind turbulence is demonstrated with calculations for a horizontal propagation 
path in the atmospheric surface layer. 

Key words: Lidar, Remote Sensing. 


1. Introduction 

Coherent Doppler lidar has become a promising new instrument for observations of the 

atmosphere. The first systems were based on the C0 2 laser 1-4 . More recently, solid state 

lasers have been successful integrated in Doppler lidar designs 5-7 . The targets for Dopper 

lidar are small atmospheric aerosol particles. Doppler lidar data can be the return from a 

single pulse or a collection of many pulses along the same line of sight. The advantages 

of multiple pulse has been discussed by Rye and Hardesty 1 * * * * * * 8 , and Frehlich and Yadlowski 9 * * * * * . 

Coherent Doppler lidar data is used to estimate the radial velocity of the targets using 

various mean-frequency estimators 8-16 . 

The actual measurement quantity is related to the statistics of the random wind field over 

the sensing volume of the pulse. The effects of the pulse volume averaging for Doppler radar 

ha s been investigated by Doviak and Zrnic 16 . The estimated Doppler frequency is defined 

in terms of the spectrum of the signal. Doppler radar data is generated by transmitting 

many pulses and sampling the returns at fixed sampling interval Ts[sec]^ which corresponds 
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to sampling the same spatial volume defined by the pulse at equally spaced time intervals. 
For this case, the statistics of the data are determined by the time evolution of the sensing 
volume of the pulse, and the spectrum is a well defined quantity. For coherent Doppler lidar, 
the data is generated from a single pulse as it travels through the atmosphere. The sensing 
volume of the Doppler lidar data is determined by the distance the pulse travels in the finite 
observation time of each estimate T = A/ 7s [sec], where M is the total number of data 
points per estimate. The spectrum is not well defined and the Doppler frequency is related 
to the statistics of the wind field over the sensing volume of the single pulse. The signal 
parameters are then related to the wind fields and lidar system parameters in a different 
fashion than for the Doppler radar. This connection will be derived to clarify the meaning 
of coherent Doppler lidar measurements and pave the way for more reliable estimators of 
the signal parameters. 

The statistical description of coherent Doppler lidar data from aerosol targets has been 
shown to be well represented as a zero mean Gaussian random process. It is a common 
practice to follow the tradition of Doppler radar and convert the real signal to the complex 
analytical signal using a quadrature receiver 16-18 or digital signal processing 19 ’ 20 . This gen- 
erates complex data with desirable statistical properties that simplify the calculation of the 
Cramer- Rao Bound 15,16,21 (the ideal performance of estimators of the signal parameters) and 
the development of efficient estimators of the signal parameters. Since the data is a Gaussian 
random process, the signal covariance provides a complete statistical description of the data. 
The signal covariance is derived for a general coherent Doppler lidar including random wind 
fields and wind shear. The derivation follows the analyses of Frehlich and Kavaya 22 for the 
performance of coherent laser radar for general atmospheric refractive turbulence. 


2. Coherent Doppler Lidar Signal 


A coherent Doppler lidar can be represented by the transmitter aperture function Wr(u), 
the transmitter laser field incident on the transmitter aperture -£<l(u,0, f)[(VFm ) ^ ]» the 
receiver aperture function VF/j(v), and the Local Oscillator (LO) field incident on the receiver 
aperture E LO {y, (^[(IFm" 2 ) 1 / 2 ]. Here, t[sec] denotes time and u[m] and v[m] denote the 
transverse two-dimensional vectors of the transmitter and receiver plane, respectively. The 
field of the transmitter laser pulse with frequency h[Hz] is written as 


E l ( u, 0, t) = U L A L (t)e L {u, 0) exp(27rt/ L f ) (1) 

where U L [J] is the laser pulse energy, ei,(u,0)[m -1 ] is the normalized spatial field which is 
assumed independent of time t (the spatial mode of the laser pulse is independent of time) 
and 


A L (t) = f f c (t')dt' 

JO 


( 2 ) 


is the normalized complex amplitude of the pulse where fc(t')[Hz] is the frequency chirp of 
the transmitted pulse and 


r \A L (t)\ 2 dt = i 

J —oo 


( 3 ) 
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The LO field in the receiver plane with frequency /lo[Hz] is given by 
Elo{v, 0 ) = Pl / Je LO {y,0)exp(2Trif LO t) 


(4) 


where Plo[W ] is the LO power, eLo( v, 0)[m l ] is the normalized spatial field which is as- 
sumed independent of time t, and 


/ OO 

|eL.Lo(u,0)|du = 1 

-OO 


(5) 


The coherent Doppler lidar typically operates under the ideal condition of quantum limited 
detection, i.e., the signal noise is dominated by the shot noise of the LO field. The effects of 
extra additive system noise is a simple extension. We will assume that the detector quantum 
efficiency ^[electrons /photon] is uniform over the surface of the detector and that all the 
LO power is collected by the detector. 

For Doppler lidar measurements of the atmosphere, the backscattered field is a zero- 
mean Gaussian random process because it is the superposition of the backscattered fields 
from many random aerosol particles. The doppler lidar signal is also a zero-mean Gaussian 
random process because it is a linear function of the backscattered field. It is convenient to 
normalize the signal by the standard deviation of the noise signal and employ the complex 
representation (Appendix A) 


Z(t) = 5(f) + N(t) = Re[z{t)\ 


( 6 ) 


where 5(0 is the normalized Doppler lidar IF current signal from the backscattered field 
and N(t) is the normalized additive detector noise. The statistical description of a real or 
complex zero-mean Gaussian random process £(<) ' s given by the covariance function 

Bt(t\,t 2 ) = H{h)((h)) (7) 

For complex random processes z(t), further simplification results if 

(z(ti)z(t 2 )) = 0 ( 8 ) 

Then, the covariance of the real random process Z(t) is related to the covariance of the 
complex random process z(t) by 


Bz{t 11^2) — 2 R e [Bz(ti, ^ 2 )] - R^\Rz{t\t ^ 2 )] 


(9) 


where Re[z] denotes the real part of the complex variable z. The function R z {ti,t 2 ) - 
B z (ti,t 2 )/2 is commonly used 9,11 15-18,21 for analysis of Gaussian random processes because 

of mathematical convenience. 

For Doppler lidar, the signal S{t) is statistically independent of the noise N(t) and 


Bz(ti,h) — Bs(ti,t 2 ) + Bx(ti,t 2 ) 


( 10 ) 


The noise N(t) is usually well described as a white-noise process over the receiver bandwidth 
B[Hz]. For the normalized Doppler lidar signal (Appendix A), 
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( 11 ) 


= R.(t,t) = SNR(t) 


where SNR is the coherent lidar signal- to- noise ratio with respect to the appropriate en- 
semble average. The Doppler lidar signal is typically uncorrelated from shot to shot because 
the aerosol targets move more than a wavelength between shots, which results in random 
phases for the backscattered signal from shot to shot. The ensemble average over the random 
phases of the backscattered fields from the aerosol particles produces (Appendix A) 


R»(hyt2) 


IqUl 

huB 


fOO roo z z 

exp{2jrt'((i - h)[&f - 2v,(p,z)/X\}dpdz 


( 12 ) 


where h — 6.626x10 4 [</s] is Planck s constant, u[Hz ] is the optical frequency, c[m/sec] is 
the^speed of light in a homogenous atmosphere, K{ p, z) is the one-way irradiance extinction, 
/?(P> 2 )[m _1 sr _I is the aerosol backscatter coefficient, t/ r (p, z) is the radial component of the 
velocity of the aerosol particles at coordinate (p, z), A/ = is the frequency difference 

between the transmitted field and the LO field, and the lidar system function <2(p,z)[m _1 ] 
is given by 


/ oo roc 

0)e£ o (v, 0)G(p; u, z)G( p; v, z)dudv (13) 

where G( p; u, z)[m ] is the Green s function for propagating the field from the transmitter 
coordinate (u,0) to the aerosol coordinate (p, zj. For Doppler lidar, the complex signal z(t) 
also satisfies Eq. (8) (see Appendix A) which facilitates data analysis and the development 
of estimation algorithms. The functions R z (t u t 2 ) and B z {t x ,t 2 ) are hermitian function, i.e., 

Rz(ti,t2) = R*(t 2 ,ti) (14) 

With a high pulse rate Doppler lidar, many shots can be transmitted while the the ran- 
dom Green s function G( p ; u, z ) (which is usually determined by the atmospheric refractive 
turbulence) and the random radial velocity fluctuations Ur(p> z) are essentially constant. For 
example, a surface layer lidar can transmit 10,000 shot per second before the signals from 
adjacent pulses interfere. Doppler lidar has the potential to investigate the statistics of wind 
fluctuations and refractive turbulence. 

The ensemble average over refractive turbulence produce 

tlnUr f°° f°° 

Rs(ti,t 2 ) = -j—g J o J_ oo A L( t I - 2z/c)A m L (t 2 - 2z/c)0(p,z) 

exp - t 2 )[Af - 2v r (p, *)/A]} dpdz (15) 

where 0( p, 2 r)[m -3 ] is the gain density 

0(p,z) = K 2 (p,z)0(p,z)c(p,z) , ( 1 6 ) 

and c(p,z)[m~ 2 ] is the coherent responsivity density 22 of the coherent Doppler lidar, i.e., 
c(p,z) = \Q(p,z )\ 2 . n 7 i 
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For many applications, the effects of refractive turbulence are negligible and 

G(p;u,*) = G J {p\u,z) = ^exp -^(p - u) 2 (18) 

where G^(p; u, z) is the free space Green’s function and k = 27r/A[m -1 ] is the wavenumber 
of the laser field. The fluctuations of the radial component of the wind field v r (p,z), the 
atmospheric extinction K(p,z), and the backscatter coefficient /?(p, z), are the remaining 
random quantities. 

Churnside and Yura 23 considered the signal covariance as an ensemble average over 
v r (p, z) by assuming v T {p,z) was a Gaussian random variable with mean value v z and 
variance <r 2 . Then, Eq. (15) becomes 

R,(ti,h) = - 2 z/e)Al(t 2 - 2z/c)0(p,z) 

exp [2jri($i — t 2 )(Af — 2v z /\) — 87r 2 <7 2 /A 2 ] dpdz (19) 

For typical lidar experiments, the analysis of Doppler lidar data is restricted to time intervals 
much smaller than the temporal scale of u r (p, z). In many cases, parameter estimation is 
performed using data from a single shot and an ensemble average over the random wind 
v r (p, z) is not appropriate. However, if the sensing volume of the pulse is much larger than 
the spatial scale of the fluctuations of v r (p, 2 ), the signal covariance has the form of Eq. 
(19). [See Fig. 6]. 

A typical Doppler lidar has a narrow beam compared to the dimensions of the sensing 
volume in the transmit direction z. Then, the function 0(p, z) behaves like a delta function 
in Eq. (15) and 

) = ^J"H(z)A L (h - 2 z/c)Al(t 2 - 2 z/c) 

exp {2 ni(ti — t 2 )[Af — 2u r (0, z)/\]} dz (20) 

where the system gain /f(z)[m -1 ] is 

H(z) = KH0,z)P(0,z)C(z) ( 21 ) 

and 

C(z)= f c(p,z)dp (22) 

J — OO 

is the coherent responsivity 22 of the coherent Doppler lidar. 

The center of a range gate for Doppler lidar measurements is defined by 

20 = cto/2 (23) 

and the total observation interval corresponding to this range gate is (to — T/2, to + T{ 2). 
If fl(p,z), K(p,z), v r (p, 2 ), and c(p, z) are constant over the sensing volume of this range 
gate, we have for all (f 1 , t 2 ) in the observation interval 
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roo 

R.(t u t 2 ) = R.(t) = SNR(zo) exp(2jrir/ m ) / A L (t)A m L (t - r)dt (24) 

where 

- !££*<«•> (25) 

is the SNR for a pulse centered at range z 0 , the mean frequency f m [Hz) of the Doppler lidar 
signal is defined as 

f m = Af-2v r (0,z o )/\ (26) 

r = tj - t 2 [sec], and the limits of integration have been extended to infinity. The signals z(t) 
and Z(t) are stationary Gaussian processes over the observation interval (t 0 — Tj 2, t 0 + T/ 2), 
i.e., the joint probability density function is only a function of the time difference r . The 
signal spectrum for a stationary real signal Z(t) is defined as 

/ CO 

Bz{r)exp(-2niTf)dT (27) 

•OO 

or 

*z(/) = SNR(zo)Sl(f — f m ) (28) 


where 

f2(/) = | f AL(t)exp{-2nift)dt\ 2 (29) 

is the spectrum of the transmitted pulse. For a stationary complex signal z(t), the spectrum 
is defined the same way as for a real signal with 5z(r) replaced with R z (r ) [see Eq. (9)]. 
The signal spectrum for aerosol targets is also the same as the signal spectrum from a rigid 
or diffuse hard target. 

For a Gaussian pulse with 1/e intensity radius <r[sec] and linear frequency chirp <f>[sec ~ 2 ] 

M‘) = ex P (" h + (30) 

and 


R s (t) = SNR(z 0 )exp(2irif m - t 2 /t p 2 ) 

where the correlation time of the lidar signal Tp[sec\ is given by 
J_ _ _j_ _L_ 

Tp 2 T A 2 + T C 2 

where r^fsec] is the correlation time due to the amplitude profile of the pulse, i.e., 
ta = 2.0 

and Tcf^ec] is the correlation time due to the frequency chirp, i.e., 

1 


(31) 

(32) 

(33) 

(34) 
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3. Effects of Wind Shear 


Doppler lidar data permits estimation of the radial component of the wind field ty(p, z) 
with a single shot. For short pulse duration, the random fluctuations of v r (p, z) can be 
approximated by the first two terms of it’s Taylor series expansion in 2 about the center of 
the range gate z 0 , i.e., 


v r (p,z) = *v(p, Zo) + g(p, Zo)(z - Zo) (35) 

where g( p, z 0 )[sec _1 ] is the wind shear of the radial component of the velocity at (p, z 0 ), i.e., 

dv r (p,z) 


9(P,z) = 


dz 


(36) 


The range dependent system gain 0( p, z ) is approximated by the first two terms of it’s 
Taylor series expansion in z, i.e., 


0(p,z) = O(p,z 0 )[l + b(p,z 0 )(z-z 0 )] 

where 6(p, z)[m -1 ] is the normalized gradient of the system gain, i.e., 
i dO(p,z) 

b(P,z = -pzrz — r — ~ 

0(p,z) dz 

When these approximations are valid and for a Gaussian pulse defined by Eq. (30) 
T ) = L0iP, *>)[! + r P b( p, 2o)(^-^ - «^-)] 


exp 


. .. , 47 ri . „.(u — fn) r r 2 

2mAfr — u r (p, z 0 )t - 2i — - — -- — 


A a Tws tt 

where n[sec\ is the centroid of and t 2 , 


dp 


P = + ^)/2, 

_ A 

TWS 2'Kr P g(p,z Q ) 

is the correlation time due to wind shear, 


(37) 


(38) 


(39) 


(40) 

(41) 


rp = caj 2 

is a measure of the range resolution of the pulse, 

1 _ 1 1 
tt 2 n 2 + t ^ 2 


(42) 


(43) 


1 11 

T\ T C T WS 


(44) 
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Note that the signal covariance depends on the centroid variable fi and is therefore non- 
stationary. For a Doppler lidar with a narrow beam compared to the dimensions of the 
range gate [see Eq. (20)] 

R a (fi,r) = SNR(z 0 )[l + 7 ( — - -«f )]exp 

a T\ 

where 

7 = 6(0, z 0 )r P = [0(0, z 0 + r P ) — 0(0, zo)]/ 0(0, zo) (46) 

is the fractional change of the gain 0(0, z) over the range resolution r P at z 0 , if the linear 
approximation of Eqs. (35) and (37) are valid. 

The effects of wind shear become important when the correlation time due to wind sheax 
Tws is less than the correlation time due to the pulse amplitude profile ta, i.e., 

ta/tws = 2irccr 2 g(0,zo)/\ > 1 (47) 

ta/tws = 0.188 for a wind shear of 0.001 sec -1 (10 m/sec per Km), A = 10 pm, and 
a = 1 fisec. 

The effects of wind shear with constant SNR [constant H(z )] over the observation interval 
are shown in Fig. 1. The wind shear produces a rapid decorrelation of the lidar signal, which 
appears as a damped oscillating signal as the centroid of the observation points /i deviate 
from the center of the observation interval to- This reflects the change in the Doppler shift 
of the signal from aerosol particles with a different velocity than the instantaneous velocity 
v r (0, zo) at the center of the range gate. The effects of a gradient in the SNR [H{z)\ is shown 
in Fig. 2. The gradient produces an amplitude change in addition to the effects of wind 
shear. 


2 TTtf m T - 2 i 


■if*- to) T 


Tws T T 


(45) 


4. Effects of Wind Turbulence 


Atmospheric wind fields are characterized by random fluctuations in space and time. For 

many cases, the statistics of the wind fields are well characterized by universal functions of 

basic parameters. For Doppler radar the effects of the random fluctuations of the wind field 

over the sensing volume has been discussed by Doviak and Zrnic 4 * * * * * * * * * * * 16 * * . For Doppler lidar, the 

situation is complicated because estimates of the velocity are computed for a single shot or 

a few shots. During the total observation period, the velocity changes little compared to the 

temporal scale of velocity flucutations. Ensemble averages of the Doppler lidar signal over 

the atmospehric velocity fluctuations were considered by Churnside and Yura 23 . The effects 

of spatial variations of the backscatter coefficient over the sensing volume of the pulse were 

investigated by Rye 24 using simulations. 

The covariance of the Doppler lidar signal for a single shot is given by Eq. (15), which 

can be written as 


R ,(tx,t 2 ) = ^ jH -2z/c)Al(U -2z/c)0(p,z) 

exp \2irif m T - ^Au r (p, z)r| dpdz 


(48) 
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where 


Aiv(p,z,z 0 ) = v r (p,z) - v r (0,z o ) 


(49) 


is the difference of the radial component of the velocity compared to v r (0, z 0 ). 

The mean frequency f m is a random variable because the wind field is random. The 
estimation error for estimates of the mean frequency f m using data from a single shot is 
difficult to predict because the signal statistics depend on the random functions Atv(p, z, zq ) 
and 0( p, z). It may be possible to estimate the random functions Av r (p, z, z 0 ) and 0( p, z ) 
as well as the other signal parameters if many shots of data can be collected during the 
temporal scale of these fluctuations. If this approach is not possible, then an ensemble 
average over the random functions Au r (p, z , zo ) and 0(p, z) is required to develop and 
analyze estimation algorithms for the signal parameters 17 . 

For well behaved atmospheric conditions, the velocity difference Av r (p, z, Zq) is approx- 
imately a zero mean Gaussian random variable. Performing the ensemble average over this 
random variable [< exp(i£) = exp(— < f 2 > /2) for a zero-mean Gaussian random variable 
£] produces 


t 2 ) = i - 2 z/c)Al(t 2 - 2z/c)0(p, z) 


exp 


2 nirfm - 


8 

A 2 


0r(P? %•) Zo) 


dpdz 


(50) 


where 


^r(p>2, Zo) — (At> r (p, 2, Z 0 ) ) 


(51) 


is the structure function of the fluctuations of the radial component of the velocity field. Note 
that the signal covariance depends on the centroid variable /i and is therefore non-stationary. 
For a narrow beam such that the function 0(p, z) dominates the p integration 


R,(tu * 2 ) = - 2z/c)Al(t 2 - 2 z/c)H(z) 


exp 2 Tirf m -p— Z) r (0, z, z 0 ) 


dz 


(52) 


For a narrow beam and a short pulse such that H(z) is constant over the sensing range of 
the pulse, 


= SNR(z 0 )J o A L (t x - 2 z/c)A m L (t 2 - 2 z/c) 


exp 


2 xir/m - 


87TV 

A 2 


£>r(0,Z,2 O ) 


dz 


(53) 


When the statistics of the velocity fields over the sensing volume of the pulse are locally 
stationary and isotropic 25 


D r { p, z, zo) = C v e 2/3 (z 0 )r 2/3 [(z - z 0 ) 2 + 4p 2 /3]/r 2 


(54) 
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where C v « 2 is the Kolmogorov constant, e[m 2 s 3 ] is the energy dissipation rate, and 
r 2 = (z — zq ) 2 +p 2 . For a Gaussian transmitted pulse with linear frequency chirp [Eq. (30)] 

r) = 2 exp(2,rir/ ” - rp( + *„) 

exp {-[(/x - to)/^ - c] 2 + 27rt><r[(/i - < 0 )/<7 - Ck - « -4/3 (C 2 + WI$)t 2 It e 2 } dqdC (55) 

where k 2 = £ 2 + g 2 , rp£ = z — z 0 , rpq = p, and the correlation time due to turbulence 
rp[sec] is given by 

— r = 87r 2 C v (crp) 2/3 / A 2 = 87r 2 £> r (0, z 0 + r P , z 0 )/ A 2 (56) 

te 

For a coherent Doppler lidar with a Gaussian transmitted spatial profile 22 
0(r F q,z) = ^ -exphrpy/er^z)] (57) 

' Ka B\ z > 

where crg{z) is the 1/e intensity radius of the transmitted beam at range z. 

For a narrow beam, the integration over q can be performed and 

r ‘^ t) = 7^!Z p b exp 

/ C(rpC + *o) exp {-[(/x - to)/* ~ C] 2 + 2iri<i><T[{p - t 0 )/a - C] T ~ C 2/3 ^ 2 / T s 2 } (58) 

J—oo K 

For a short pulse such that C(rp£ + z 0 ) is constant over the £ integration, 

R .(v, r) = S -~~ exp ^2jrtr/ m - 

f exp{-[(p-t 0 )/(T-(} 2 + 2iri<f><7[(ti-t 0 )/<T-(]T-( 2/3 T 2 /T E 2 }d{ ( 59 ) 

J - 00 W 

The effects of wind turbulence become important when the correlation time due to tur- 
bulence t e becomes less than the correlation time due to the pulse amplitude profile r^, i.e., 
when 

ta/te = 47rC„ 1/2 <r(crp) 1/3 /A > 1 (60) 

For typical stable boundary layer turbulence, e = 0.01 m 2 /sec 3 . Then r^/rp = 2.88 for 
A = 10 pm, and a = Ipsec. 

The effects of wind turbulence and the finite transverse dimensions of the beam are 
shown in Fig. 3. The decrease in the correlation scale due to turbulence is more pronounced 
for a wide beam (erg > rp) because of the increase in the Doppler shifts from the random 
fluctuations of the wind in the transverse dimensions. 

The effects of wind turbulence for the common case of a narrow transmitted beam com- 
pared to the range resolution (erg < rp) is shown if Fig. 4. When t e < ta , the signed 
covariance is reduced due to the wind turbulence. 
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For typical surface layer measurements, the statistical description of the velocity fluc- 
tuations must include the spatial scale of the turbulent wind field (sometimes called the 
outer-scale of turbulence) that is proportional to z/f[m], the height above the surface. The 
spectra of velocity fluctuations satisfy universal scaling functions. For a horizontal propaga- 
tion direction along the mean velocity, the universal spectrum of velocity flucutations in the 
direction of the mean velocity results in a universal, locally stationary, structure function of 


velocity fluctuations. For neutral stability 26 

D r { 0, z 0 + r, z 0 ) = 2<7 2 A(r/ zp) (61) 

where a r [m/sec ] is the standard deviation of the radial component of the velocity flucutations 
and the universal function 

A(x) = (ax) 2 / 3 [l + (ai) p ]' 2 ^ 3p * (62) 

where a = 0.26278 and p = 1.1948. The function A(x) is shown in Fig. 5. For small r, Eq. 
(61) is equal to Eq. (54) with p = 0 and 

c = (2<7 2 /C' v ) 3 /2 q / zp (63) 

Substituting Eqs. (30) and (61) into Eq. (53) produces 

> SNR(z 0 ) , r ! \ 

R s(^ r) = ;rl/2 exp \2wtT f m - — j 

J exp {-[(p - t 0 )/cr - C ] 2 + 27 ri<t>cr[{p - t 0 )/ cr - C]r - A {r P (J zp)r 2 /rp 2 } d(, (64) 

where 

-i- = 167 t 2 <t 2 /A 2 (65) 

th 


There are two limiting cases of Eq. (64). When the sensing volurpe of the pulse rp is less 
than the outer scale of turbulence ( Lq « zp) the integration in Eq. (64) is dominated by the 
power law region of A(x) and Eq. (64) is equal to Eq. (59). When rp is much larger than 
zh , then A (rp^/zn) « 1 and a result similar to Eq. (19) is produced. If the mean frequency 
f m is defined using the average radial velocity v z instead of the instantaneous radial velocity 
*>r(0 ,zq), then Eq. (19) is produced when rp » zp because D r (p, z,zq) — ► <7 2 instead of 
2 <t 2 [see Eq. (61)]. 

The signal covariance Eq. (64) is shown in Fig. 6 for a typical nightime condition with 
a horizontal path at height zp = 2m, a velocity fluctuation <r r = 0.5 m/sec and operating 
wavelengths of 2 pm and 10pm. This corresponds to an energy dissipation rate of t = 
0.0164m 2 /s 3 . For these parameters, the correlation scale for turbulence Tfj assuming rp >> 
zh is 1.59psec for the 10pm lidar and 0.318psec for the 2 pm lidar. For this regime, the 
signal covariance is determined by the transmitted pulse when < 77/ . 
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5. Summary 


The signal covariance for coherent Doppler lidar data is derived for general conditions for 
both real and complex data. The covariance is defined with respect to a given ensemble 
average. The most general case is the ensemble average over the random phases of the 
backscattered fields from the aerosol particles [Eq. (12)]. The covariance parameters are 
a function of the instantaneous wind field and atmospheric conditions. In this case, the 
covariance can be estimated with multiple shots over a time interval that is less than the 
temporal scale of the fluctuations in the atmosphere. The next level of ensemble average 
is over the atmospheric refractive turbulence [see Eq. (15)]. In many cases, the effects of 
refractive turbulence are negligible and Eq. (18) is valid. Under the ideal limit of uniform 
velocity, atmospheric parameters, and system parameters over the sensing volume of the 
pulse, the Doppler lidar signal is stationary and determined by the transmitted pulse Eqs. 
(24), (28), and (29). When there is a linear wind shear and linear system gain over the sensing 
volume of the pulse, the signal covariance is non-stationary (depends on n) and given by 
Eq. (39) and for a narrow lidar beam by Eq. (45). The signal parameters are expressed in 
terms of the instantaneous atmospheric variables and lidar system parameters. The mean 
frequency f m is defined in terms of the instantaneous radial velocity in the center of the range 
gate of interest for a single shot. This permits analysis and development of more efficient 
estimation algorithms which are valid when the pulse sensing volume is small compared to 
the variations of the instantaneous atmospheric variables. The coherent lidar signal is not 
stationary over the observation interval and spectral domain estimators will be less efficient 
for this case because they assume stationary data. When the atmospheric parameters are 
unable to be modeled over the observation interval (e.g., velocity fluctuations are not a 
linear function of space over the sensing volume of the pulse), the analysis and development 
of estimation algorithms requires the ensemble average over the unknown parameters 18 . The 
signal covariance for this case with general statistics for the velocity fluctuations is given by 
Eqs. (50), (52), and (53). For isotropic and homogeneous wind fluctuations described by the 
Kolmogorov turbulence spectrum over the sensing volume of the pulse, the signal covariance 
is given by Eqs. (50), (52), and (53) with Eq. (54). and for a Gaussian pulse defined by Eq. 
(30), the correlation scale is given by te Eq. (56). When the sensing volume of the pulse 
is much larger than To, then Eqs. (50), (52), and (53) are valid with D r (p, z,zo) = 2<x r 2 
is valid and the correlation time is given by 77/ Eq. (65). For intermediate regimes, the 
universal description of the velocity structure function is required [see Eq. (61)] For the case 
of horizontal propagation in the atmospheric surface layer and typical nightime conditions, 
the correlation time due to the effects of turbulence are less than 2 fxsec for a 10 pm lidar 
and less them 0.2 ftsec for a 2 micron lidar. The correlation time will be reduced even more 
for typical daytime unstable convective conditions. For wind shear and wind turbulence, 
the correlation time is proportional to the operating wavelength A [see Eqs. (41), (56), (65)] 
and the critical pulse length where the effects of the wind fluctuations dominate over the 
contribution of the pulse is also proportional to A. 
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Appendix A: Coherent Doppler Lidar Signals 

The total coherent or heterodyne Doppler lidar signal ij(f)[A] is the sum of the Doppler IF 
signal current ts(f)[A] an d the detector noise current t/v(f)[/4] which is conveniently written 
in complex notation as ir(i) = f?e[j/(f)], i^(t) = f?e[u(f)], where 

y(t) = u(t ) + J^t] Q (w)E s (w, L,t)E m LO (w, L)exp(2niAft + i0 s )dv> (Al) 

Go is the dimensionless amplifier gain, e = 1.602 x 10~ l9 [C / electron] is the electronic 
charge, rjq(w)[electrons / photon] is the detector quantum efficiency function on the detector 
surface, £'s(w, L, t ) is the backscattered field on the detector surface, Elo{w, L) is the LO 
field on the detector surface, 9s[rad\ is the random phase of the backscattered field compared 
to the LO field, Af[Hz\, is the frequency difference between the transmitted and LO field, 
w[m] is the transverse coordinate on the detector surface and <fw denotes two-dimensional 
integration over the detector surface. For ideal optical systems and for a infinite uniform 
detector, the fields on the detector can be represented in terms of the fields in the plane of 
the receiver optics, i.e., 

/ oo 

Es{v, 0, t)E* LO (\, 0) Wfl(v) exp(2riAft + i0 s )dv (A2) 

-OO 

where Es(v, 0, t) and £lo(v, 0) are the backscattered field and LO field in the receiver plane. 

If the signal noise is dominated by the shot noise generated by the local oscillator field 
(quantum limited detection), the average noise power is 

i 2 N (t) = i«(tK(0 = 2G D 2 e 2 Br }Q P L ol(hu ) (A3) 

where B[Hz] is the noise bandwidth. It is convenient to normalize the signals by the rms of 
the noise current, i.e., 

z(t) = y(t)/i 2 N (t)' /2 = s{t) + n(t). (A4) 

The backscattered field from a point scatterer at coordinate (p, z) with scattering cross 
section crs is 

Bs(v,0,t) = \(t)/ 2 K(p,z) [ E L (u,0,t - 2z/c)W T (u)G(p;\i,z) 

G(p; v, z) exp[i0(p, z ) — 43rtfu r (p, z)/A]du (A5) 

where 0(p,z) is the random phase of the backscattered field. The total backscattered 
field is the sum of all the backscattered fields from the aerosol particles illuminated by 
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the transmitted field. Substituting the total backscattered field into Eqs. (A2), (A4), (9), 
(1), (4), performing the ensemble average over the random phases 0(p, z) using the iden- 
tity exp[i(0(p, zi) + 0(p, z-i)] = 0, and converting the summation over aerosol particles to 
integration 22 produces Eq. (12). The same procedure produces Eq. (8). 
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Figure 1 

Fig. 1. Effects of wind shear on the normalized signal covariance [Eqs. (9), (30), (45), 
■y — f m = <^=0] as a function of t/t a with constant SNR over the observation interval and 

( ) n = t 0 , no wind shear, (....) fx-t 0 = cr and t a = tws\ ( ) fi-t 0 = 2(r and t a = nvs; (- • 

- . ) fi — to — <7 and t a — 2 rws! an d (••- •• -) M — = and t a = 2 t\vs- 





B s ( n , x )/B s ( t 0 , 0 ) 


Figure 2 


Fig. 2. Effects of wind shear and gradients in SNR over the observation interval on the normalized 
signal covariance [Eqs. (9), (30), (45), f m = <j>= 0] as a function of t/t a compared with constant 

SNR over the observation interval and no wind shear ( ). The curves represent (....) n - 1 0 = a, 

ta = nvs* and 7 = 0.2; ( ) /x -t 0 = cr, t a = rws, and 7 = 0.4; (-.-.) fi — t 0 = 2 <r, t a — 2 Tws, 

and 7 = 0; (.. - -) n - t 0 = 2 cr, t a = 2 tws, and 7 = 0.2; and (..- - -) fi-t Q = 2(7, r A = 2 tws, 

and 7 = 0.4. 




B s ( |A , T )/Bg ( t„ , 0 ) Bg ( (i. , X )/Bg ( t 0 , 0 ) 


riyure 3 


Fig. 3. Effects of the finite transverse dimensions of a lidar beam and wind turbulence on the 
normalized signal covariance [Eqs. (9), (30), (55), f m = 4>=0] as a function of r/r A with constant 
SNR over the observation interval and case a) fi = t 0 with ( — ) no wind turbulence; (...) t e - t a / 2, 
and o B « rp; (- -) t e = t a /2 , and o B = r P / 2; (•-•-) T E = t a /2, and o B = r P ; (- - . - - .) 
r E = t a /2, and o B = 2 r P \ and case b) the same as case a) but with n - t 0 = a. 



t/t a 



( o ‘ °» ) s a/( i ‘ ) s a 


Figure 4 


Fig. 4. Effects of wind turbulence on the normalized signal covariance [Eqs. (9), (30), 
(59), f m = 4 >= 0 ] as a function of t/t a with constant SNR over the observation interval and 
narrow transmitted beam compared to the range resolution ( erg « r p) for case a) /j = <o with 
( — ) no wind turbulence; (...) te = t a /2\ (- -) te = t a \ and case b) is the same as case a) but 
with fi — t 0 = cr. 







f-igure b 

Fig. 6. Normalized signal covariance [Eqs. (9), (30), (64), f m = <j >= 0 ] as a function of r for 
fi-t 0 and a 10 fim and 2/xm lidar for a horizontal propagation path at a height z„ = 2m in the 
atmospheric surface layer under neutral stability with a standard deviation of the radial component 
of velocity fluctuations a T = 0.5 m/sec. The correlation scale assuming the lidar pulse dominates 
the signal statistics r A = 2a is (_) 0.2 nsec, (...) 0.4/wee, (- -) 1.0/«ec, (.-.-) 2.0 usee, (- .. - ) 

4.0 nsec, and (- - . - - .) 10.0/xsec. ’ ’ 
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